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Pulsars are wonderful gravitational probes. Their tiny size and stellar mass give their 
rotation periods a stability comparable to that of atomic frequency standards. This is es¬ 
pecially true of the rapidly rotating “millisecond pulsars” (MSPs). Many of these rapidly 
rotating pulsars are in orbit with another star, allowing pulsar timing to probe relativis¬ 
tic perturbations to the orbital motion. Pulsars have provided the most stringent tests of 
theories of relativistic gravitation, especially in the strong-field regime, and have shown 
that Einstein’s general theory of relativity is an accurate description of the observed 
motions. Many other gravitational theories are effectively ruled out or at least severely 
constrained by these results. MSPs can also be used to form a “Pulsar Timing Array” 
(PTA). PTAs are Galactic-scale interferometers that have the potential to directly de¬ 
tect nanohertz gravitational waves from astrophysical sources. Orbiting super-massive 
black holes in the cores of distant galaxies are the sources most likely to be detectable. 
Although no evidence for gravitational waves has yet been found in PTA data sets, the 
latest limits are seriously constraining current ideas on galaxy and black-hole evolution 
in the early Universe. 
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1. Introduction 

Pulsars are rotating neutron stars that emit beams of radiation which sweep across 
the sky as the star rotates. A beam sweeping across the Earth may be detected as a 
pulse that repeats with a periodicity equal to the rotation period of the star. Because 
of the large mass of neutron stars, ~ 1.4 M©, and their tiny size, radii ~ 15 km, (see 
Ref. [76] for a review of neutron-star properties) the rotation period of neutron stars 
is incredibly stable, with a stability comparable to that of the best atomic clocks 
on Earth. This great period stability, combined with the fact that pulsars are often 
in a binary orbit with another star, makes them wonderful probes of relativistic 
gravity. Tiny perturbations to their period resulting from, for example, relativistic 
effects in a binary orbit, may be detected and compared with the predictions of a 
gravitational theory. Pulsars may also be used as detectors for gravitational waves 
passing through the Galaxy. To separate the effects of gravitational waves from 
other perturbations, signals from pulsars in different directions on the sky must be 
compared - exactly analogous to the way laser-interferometer gravitational-wave 
detectors compare laser phases in perpendicular arms. 


1 






April 21, 2015 0:29 WSPC/INSTRUCTION FILE grlOO'rnm'ijmpd 


2 


More than 2400 pulsars are now known If] The vast majority of them reside 
in our Galaxy, typically at distances of a few thousand light-years from the Sun. 
Their pulse periods, P, range between 1.4 milliseconds and 12 seconds and fall into 
two main groups. The so-called “normal” pulsars have periods longer than about 
30 milliseconds and the “millisecond” pulsars (MSPs) have shorter periods. MSPs 
comprise about 15% of the known pulsar population. 

Although pulsar periods are very stable, they are not constant. In their own 
reference frame, all pulsars are slowing down, albeit very slowly. Pulsars are pow¬ 
ered by their rotational energy. They have extremely strong magnetic fields and, 
as they spin, they emit streams of relativistic particles and so-called “dipole radi¬ 
ation” , electromagnetic waves with period equal to the rotation period of the star. 
These carry energy away from the star resulting in a steady increase in the pulse 
period. Typical rates of period increase, P, are a part in 10 15 for normal pulsars and 
much less for MSPs. Assuming that the surrounding magnetic field is predominantly 
dipole, the characteristic age r c of a pulsar is given by 

Tc = P/ (2P) (1) 

and their surface dipole magnetic field strength (in gaus:JH) is 

B s « 3.2 x 10 19 (PP) 1/2 G. (2) 

Figure [T] shows the distribution of pulsars on the P - P plane, with several different 
types of pulsars indicated. For most normal pulsars, r c is between 10 3 and 10 7 years 
and B s is between 10 11 and 10 13 G. For MSPs, the corresponding ranges are 10 9 to 
10 11 yrs and 10 s to 10 10 G. 

Normal pulsars and MSPs have quite different evolutionary histories. Most if not 
all normal pulsars are formed in supernova explosions at the death of a massive star. 
They age with relatively constant B s until the pulse emission mechanism begins to 
fail when r c reaches about 10 6 yrs. Many young pulsars (r c < 10 4 yrs) are located 
within supernova remnants, with the most famous example being the Crab pulsar, 
PSR B0531+21, located near the centre of the Crab Nebula. Most young pulsars lie 
relatively close to the Galactic Plane, consistent with the idea that they are formed 
from massive stars. MSPs on the other hand are much more widely distributed in the 
Galactic halo. They are believed to originate from old, slowly rotating and probably 
dead neutron stars that accrete matter and angular momentum from an evolving 
binary companion. This “recycling” process increases their spin rate so that they 
have periods in the millisecond range and re-energises the beamed emission (see, 
e.g., Ref. fTTj) . 

Figure |T] shows that the majority of MSPs are members of a binary system, 
consistent with this formation scenario. The accretion also suppresses their apparent 
magnetic field, so that MSPs have P values about five orders of magnitude less than 


a See the ATNF Pulsar Catalogue: http://www.atnf.csiro.au/research/pulsar/psrcat 
b l gauss (G) = 10 —4 tesla 
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normal pulsars. Since the level of intrinsic period irregularities is related to /hLLLZJ 
it is this low B s that makes MSPs extremely stable clocks and suitable tools for 
the study of relativistic gravitation. Figure [Q also indicates the class of radio-quiet 
pulsars. Essentially all of the 72 known radio-quiet pulsars are young and solitary. 
They include the so-called “magnetars” which lie in the upper right side of the P - 
P diagram where dipole magnetic fields are strongest. 

Most double-neutron-star systems lie in the zone between the MSPs and the 
normal pulsars. These systems are believed to have been partially recycled prior 
to the formation of the second-born neutron star. In almost every case, the ob- 



Fig. 1. Distribution of pulsars on the P — P plane with radio-loud and radio-quiet pulsars indi¬ 
cated. Binary pulsars are indicated by a circle around the point and, for those with a neutron-star 
companion, the circle is filled in. Lines of constant characteristic age (r c ) and surface-dipole mag¬ 
netic field ( B s ) are shown. 
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served pulsar is the recycled one since it has a much longer active lifetime than the 
newly formed young pulsar. A famous exception to this is the Double Pulsar (PSR 
J0737—3039A/B) in which the second-born star (B) is (or was) still visible.^ On 
Figure [T] the B pulsar is the solitary double-neutron-star system on the right side 
of the plot. The double-neutron-star system identified near the middle of the plot 
contains the relatively young pulsar, PSR J1906+0746. There is some doubt about 
the nature of the companion in this system it could be a heavy white dwarfd 39 ^ 
These double-neutron-star systems and their use as probes of relativistic gravity 
are discussed in some detail in Section [2711 below. 

1.1. Pulsar Timing 

The most important contributions of pulsars to investigations of gravitational the¬ 
ories and gravitational waves rely on precision pulsar timing observations. These 
allow both the relativistic perturbations to binary orbits to be studied in detail and 
the potential detection of the tiny period fluctuations generated by gravitational 
waves passing through our Galaxy. Because of the importance of pulsar timing to 
these studies, we give here a brief review of its basic principles. 

The basic observable in pulsar timing is the time of arrival (ToA) of a pulse at an 
observatory. In fact, because of signal/noise limitations and the intrinsic fluctuation 
of individual pulse shapes, ToA measurements are based on mean pulse profiles 
formed by synchronously averaging the data, typically for times between several 
minutes and an hour. The time at which a fiducial pulse phase (usually near the pulse 
peak) arrives at the telescope is determined by cross-correlating the observed mean 
pulse profile with a standard pulse template. A series of these To As is measured 
over many days, months, years and even decades for the pulsar of interest. 

These observatory ToAs are affected by the rotational and orbital motion of 
the Earth (and for satellite observatories, the orbital motion of the satellite). To 
remove these effects, the observed To As are referred to the solar-system barycentre 
(centre of mass) which is assumed to be inertial (unaccelerated) with respect to 
the distant Universe])] This correction makes use of a solar-system ephemeris giving 
predictions of the position of the centre of the Earth with respect to the solar- 
systern barycentre. Such ephemerides, for example, the Jet Propulsion Laboratory 
ephemeris DE 42 1,03 are generated by fitting relativistic models to planetary and 
spacecraft data. The correction also takes into account the relativistic variations in 
terrestrial time resulting from the Earth’s motion. 

The resulting barycentric ToAs are then compared with predicted pulse ToAs 
based on a model for the pulsar. The pulsar model can have 20 or more parameters; 
generally included are the pulse frequency (v = 1/P), frequency time derivative 
(!'), the pulsar position and the five Keplerian binary parameters if the pulsar is a 


c This neglects any acceleration of the solar-system barycentre resulting from, for example, Galactic 
rotation. For some precision timing experiments, such effects are taken into account at a later stage 
of the analysis. 
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member of a binary system. If the data are recorded at different radio frequencies, 
it is also necessary to include the dispersion measure (DM) which quantifies the 
frequency-dependent delay suffered by the pulses as they propagate through the 
interstellar medium. 

The differences between the observed and predicted ToAs are known as “timing 
residuals”. Errors in any of the model parameters result in systematic variations in 
the timing residuals as a function of time. For example, if the model pulse frequency 
is too small, the residuals will grow linearly as illustrated in Figure [2] The required 
correction to the pulse frequency is given by the slope of this variation. An error in 
the pulsar position results in an annual sine curve which arises from the barycentric 
correction. The phase and amplitude of this curve give the corrections to the two 
position coordinates. Similarly, a pulsar proper motion results in a linearly growing 
sine curve (away from the reference epoch). For a sufficiently close pulsar, the cur¬ 
vature of the wavefront results in a biannual term in the residuals, and offsets in 
binary parameters result in terms which vary with the orbital periodicity. 

Since offsets in different parameters in general result in different systematic 
residual variations, it is possible to do a least-squares fit to the observed residuals 
for the correction to any desired parameter^] The precision of the fitted parameters 
is often amazingly high. For example, the relative precision of the pulse frequency 
determination is ~ St/T , where St is the typical uncertainty in the ToAs and T is 
the data span. Since, for MSPs at least, St is often < 1 ys and T is often many years, 
the relative precision of the measured v can easily exceed 1:10 14 . Similarly pulsar 
positions can be measured to micro-arcseconds and binary eccentricities measured 
to 1 : 10 8 . These very high precisions often allow higher-order terms, for example, 
resulting from relativistic perturbations to the binary orbit, to be measured. 

2. Tests of Relativistic Gravity 

2.1. Tests of General Relativity with Double-Neutron-Star Systems 

2.1.1. The Hulse-Taylor binary, PSR B1913+16 

The discovery at Arecibo Observatory in 1974 of the first-known binary pulsar, 
PSR B1913+16, by Hulse and Tayloi® was remarkable in a number of respects. 
Firstly, it showed that pulsars with short pulse periods but large characteristic ages 
(59 ms and 10 8 yr, respectively, for PSR B1913+16) could exist. The period of 
PSR B1913+16 was second only to the Crab pulsar, but its age was enormously 
greater than that of other short-period pulsars known at the time. Secondly, it was 
in a binary orbit with a relatively massive star, very likely another neutron star, 
showing that an evolutionary pathway to such systems existed. Thirdly, its orbital 
period was extraordinarily short, only 7.75 h, its eccentricity large, ~ 0.617, and, 
as shown in Figure [3l its maximum orbital velocity very high, ~ 300 km s _1 or 


d The data span must be sufficiently long to avoid excessive covariance between the variations for 
different fitted parameters. 
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Fig. 2. Variations in timing residuals for offsets in several parameters. The model pulsar has a 
pulse frequency of 300 Hz (3.3 ms period) and is in an eccentric (e = 0.4) binary orbit of period 
190 days. The reference epochs for period, position and binary phase are at the middle of the 
plotted range. 


0.1% of the velocity of light. As was immediately recognised by Hulse and Taylor, 
these properties opened up the possibility that relativistic perturbations to the orbit 
were potentially measurable. Lowest-order relativistic effects go as (v/c) 2 , and so the 
variations are of order 1:10 6 , enormous by the standards of pulsar measurements. 

Relativistic effects in binary motion can be expressed in terms of “post- 
Keplerian” parameters that describe departures from Keplerian motion (see, e.g., 
Ref. I129[) . The first such parameter to be observed was periastron precession. 135 
In Einstein’s general theory of relativity (GR) the rate of periastron precession 
(averaged over the binary orbit) is given by: 

w = 3(P 6 /27 t)- 5 / 3 (T 0 M) 2 / 3 (1 - e 2 )" 1 (3) 

where ui is the longitude of periastron measured from the ascending node, Pb is the 
orbital period, T 0 = GMq/c 3 = 4.925490947/zs, G is the Gravitational Constant, 
Mq is the mass of the Sun, e is the orbital eccentricity and M = m\ + m 2 , the 
sum of the pulsar mass m\ and the companion mass m 2 in solar units. It is worth 
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Fig. 3. Velocity curve for the Hulse-Taylor binary pulsar, PSR B1913+16. (Ref. Ifilt 


noting that this is the same effect as the excess perihelion advance of Mercury that 
was used by Einstein‘S as an observational verification of GR. The relativistic effect 
for Mercury is just 43 arcsec per century, minuscule compared to the 4?22 per year 
predicted and observed for PSR B1913+16. 

The next most significant parameter, normally labelled 7 , describes the combi¬ 
nation of gravitational redshift and 2nd-order (or transverse) Doppler shift, both of 
which have the same dependence on orbital phase. In GR 

7 = e(P h /27r) 1 / 3 TQ /3 M _4/3 m 2 (mi + 2m 2 ). (4) 


Since the Keplerian parameters are very well determined, measurement of u> and 7 
gives two equations in two unknowns, m 1 and m 2 , and so the two stellar masses can 
be determined. For PSR B1913+16, these are both close to 1.4 M 0 , confirming the 
double-neutron-star nature of the system. An important consequence of this was 
that the two stars could safely be treated as point masses in GR, thereby allowing 
precise tests of the theory. 

Given the Keplerian parameters and the two masses, the next post-Keplerian 
parameter, orbital decay due to the emission of gravitational waves from the system, 
given in GR by 


i\ = - 


192tt 

5 



, 73 2 37 . 

1 4-e 2 H-e 

24 96 


(1 - e 2 ) -7/2 T | /3 mim 2 M- 1 / 3 , 


(5) 


could be predicted and compared with observation. This therefore constituted a 
unique test of GR in strong-field gravity. After just four years, the Pb term was well 
measured and shown to be consistent with the GR prediction!-^ The orbital decay, 
including more recent results, is illustrated in Figurc[I] the ratio of observed to pre¬ 
dicted orbit decay is 0.997 ± 0.002^^21 It should be mentioned that this near-perfect 
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agreement relies on correcting the observed Pj, for the differential acceleration of 
the solar system and the pulsar system in the gravitational field of the Galaxy. This 
correction is —0.027 ± 0.005 x 10 ~ 12 or about 1% of the observed value, and the 
uncertainty in the final result is dominated by the uncertainty in this correction. 
Unfortunately, it is unlikely that this uncertainty will be significantly reduced in 
the near future since it mainly depends on the poorly known distance to the binary 
system. 



Fig. 4. Comparison of the observed and predicted orbital period decay in PSR B1913+16. The 
orbital decay is quantified by the shift in the time of periastron passage with respect to a non¬ 
decaying orbit. The parabolic curve is the predicted decay from GR. (Ref. [T43 1 


Two more post-Keplerian parameters, denoted by r and s, relate to the Shapiro 
delay suffered by the pulsar signal while passing through the curved spacetime 
surrounding the companion starP 23 The relations for them in GR are as follows: 

r = Tqto 2 (6) 

/ p \ —2/3 

s = sin/ = (ai/c) sin/ ( — J Tq ' 3 (7) 

where (ai/c) sin / is the projected semi-major axis of the pulsar orbit (in time units), 
a Keplerian parameter. The Shapiro delay is generally only observable when the 
orbital inclination is relatively close to 90°, that is, the orbit is seen close to edge- 
on. For the Hulse-Taylor binary pulsar the orbital inclination is about 47°and the 
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Shapiro delay is small and covariant with Keplerian parameters. However it has 
been observed in a number of other neutron-star binary systems as will be discussed 
below. 

Another relativistic effect, geodetic precession, has observable consequences for 
PSR B1913+16 and several other binary pulsars. A neutron star formed in a super¬ 
nova explosion receives a “kick” during or shortly after the explosion which typically 
gives the star a velocity of several hundred km s _1 .^If the pulsar is a member of a 
binary system which is not disrupted by the kick, the orbital axis is changed so that 
pulsar spin axis, which was most likely aligned with the orbital axis prior to the 
explosion, is no longer aligned. Since kick velocities are often comparable to or even 
larger than the orbital velocities, the misalignment angle can be large. Precession 
of the spin axis will therefore alter the aspect of the radio beam seen by an observer 
and may even move the beam out of the line of sight. 

In GR, the precessional angular frequency is given by 


i — I f 5/3 r 2/3 f»2(4mi + 3m 2 ) 
p ~ 2 \2n J 0 (1 — e 2 )M 4 / 3 


( 8 ) 


For PSR B1913+16, the corresponding precessional period is about 297 years, so the 
aspect changes over observational data spans are small. Never-the-less changes in 
the relative amplitude of the two peaks in the PSR B1913+16 profile were reported 
by Weisberg et alP^H and attributed to the effects of geodetic precession with a 
“patchy” beam. For a basically conal beam geometry, the separation of the two 
profile components would be expected to change and evidence for this was first 
found by Kramei ^ 2 leading to an estimate of the misalignment angle of about 22°. 
A data set extending to 2001 was analysed by Weisberg and Taylorp^ obtaining 
results similar to those of Kramer . 72 Their best-fitting model gives the “peanut” 
shaped beam shown in Figure [5] Clifton and Weisber^l have shown that a set of 
circular nested emission cones can also give apparent pulse-width variations similar 
to those observed. Over the 20-year interval covered by the data set, the “impact 
parameter” (minimum angle between the beam symmetry axis and the observer’s 
line of sight) has changed by a rather small amount, about 3°, from —3?5 to —6?5. 
Extrapolation of this model suggests that the pulsar will become unobservable in 
about 2025. While this result is compatible with relativistic precession, it is not 
possible to derive an independent measure of the precessional rate from these data. 


2 .1.2. PSR B1534+12 

PSR B1534+12, discovered by Wolszczan in 1990p-^l is a binary system with 
parameters quite similar to those of B1913+16, notably a short orbital period 
(~ 10.1 hours), relatively high eccentricity (~ 0.27) and a compact orbit about 
60% larger than that of PSR B1913+16. Analysis of less than a year’s data al¬ 
ready allowed measurement of two post-Keplerian parameters, w and 7 , thereby 
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Fig. 5. Contours of the symmetric part of the radio emission beam for PSR B1913+16, obtained 
by fitting to the variation in pulse width over the interval 1981 to 2001. (Ref. fl45tl 


determining the masses of the two stars and confirming that they are both neu¬ 
tron stars. These results also showed that the orbit was more edge-on than that of 
PSR B1913+16, with an implied inclination angle of about 77°. Analysis of timing 
data extending over 22 years by Fonseca et alP^l built on earlier results by Stairs 
et al.pAL with significant detections of r and s, the Shapiro-delay parameters, and 
the orbital decay, Pb, giving five post-Keplerian parameters and three independent 
tests of GR. A fourth test of GR, albeit less precise, was provided by an analysis of 
the evolution of the profile shape and polarisation, yielding a rate for the geodetic 
precession of the pulsar spin axis fl p = 0?59(j"og 12 yr _1 which is consistent with the 
value predicted by GR. Figure [ 6 ] shows the so-called “mass-mass” diagram for PSR 
B1534+12, illustrating these constraints. 

If GR gives a correct description of the post-Keplerian parameters, all of these 
constraints should be consistent with an allowed range of m\ and nri 2 - For PSR 
B1534+12, the masses are most accurately constrained by Cj and 7 but the predicted 
constraint on Pb appears to be inconsistent. As for PSR B1913+16, the measured 
value of Pb must be corrected for kinematic effects resulting from the differential 
acceleration of the binary and solar system in the Galaxy. PSR B1534+12 is much 
closer to the Sun that PSR B1913+16 and the so-called “Shklovskii” effectan 
apparent radial acceleration resulting from transverse motion of the binary system 
(P/P = [Pd/c, where \x is the pulsar proper motion and d is the pulsar distance), 
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is also important. The distance estimate is based on the pulsar DM and a model 
for the free electron density in the Galaxy and is quite uncertain. Stairs and her 
colleague^- 133 ® inverted the problem, assuming that the GR prediction for P& is 
correct, thereby deriving an improved value for the pulsar distance, a technique 
first suggested by Bell and BailesP^ 



Pulsar mass (M 0 ) 


Fig. 6 . Plot of companion mass (m 2 ) versus pulsar mass (mi) for PSR B1534+12. Constraints 
derived from the measured post-Keplerian parameters assuming GR are shown pairs of lines with 
separation indicating the uncertainty range. For 6j and 7 , the uncertainty ranges are too small to 
see on the plot. In addition, a constraint from the measured spin precession rate is shown. (Ref.[49} 


2.1.3. The Double Pulsar, PSR J0737-3039A/B 

The discovery of the Double Pulsar systeirl 21 ^ heralded a remarkable era for in¬ 
vestigation of relativistic effects in double-neutron-star systems. In this system, the 
A pulsar was formed first and subsequently spun up to approximately its current 
period of 23 ms by accretion from its evolving binary companion. The companion 
then imploded to form the B pulsar which has since spun down to its current period 
of about 2.8 s. The orbital period is only 2.5 hours, less than a third of that for 
PSR B1913+16, and the projected semi-major axis aisin i is about 60% that of 
PSR B1913+16. These parameters imply relativistic effects much larger than those 
seen for PSR B1913+16; for example, the predicted relativistic periastron advance 
is 16?9 yr , more than four times the value for PSR B1913+16. Added to that, the 
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orbit is seen within a degree or so of edge-on, not only allowing detailed measure¬ 
ment of the Shapiro delay, but also resulting in eclipses of the A pulsar emission by 
the magnetosphere of the B pulsar. Finally, the still-unique detection of the second 
neutron star as a pulsar allows a direct measurement of the mass ratio of the two 
stars from the ratio of the two non-relativistic Roemer delays (ai sin i) Q 

Four post-Keplerian parameters (u>, 7 , r and s) were detected in just seven 
months of timing data from the Parkes 64-m and Lovell 76-rn Telescopes.®^ Further 
observations, including data from the Green Bank 100-m Telescope, with a 2.5-year 
data span give the currently most stringent test of GR in strong-field conditionsP 1 7 'I 
Figure [ 8 ] shows the mass-mass diagram based on these results together with the 
measurement of geodetic spin precession for the B pulsar described below. A total of 
six post-Keplerian parameters together with the mass ratio R gives five independent 
tests of GR. As well as the three post-Keplerian parameters, u>, 7 and Pf, observed 
for PSRs B1913+16 and PSR B1534+12 (Figure [ 6 ]), for the Double Pulsar we have 
the mass ratio R , the Shapiro delay parameters r and s and a measurement of the 
rate of geodetic precession fi p . The observed Shapiro delay (Figured]) shows that 
the J0737— 3039A/B orbit inclination angle is 88?7 ± 0?7. 



120 180 240 

Longtiude (deg) 


Fig. 7. Observed Shapiro delay as a function of orbital phase for PSR J0737—3039A. Upper panel: 
timing residuals after fitting for all parameters except the Shapiro-delay terms, r and s, with these 
set to zero. Lower panel: the full Shapiro delay obtained by taking the best-fit values from a full 
solution, but setting r and s to zero. The red line is the prediction based on GR. (Ref. [ 74)1 


As mentioned above, this nearly edge-on view of the orbit results in eclipses 
of of the A-pulsar emission by the magnetosphere of the B pulsar. These eclipses 


e The B pulsar became undetect able as a radio pulsar in 2008, most probably because the beam 
precessed out of our line of sightBecause of uncertainty about the beam shape, the date of its 
return to visibility is very uncertain, but it should be before 2035. 
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last only about 30 s, showing that the B-pulsar magnetosphere is highly modified 
by the relativistic wind from pulsar A. 8 ^ Remarkably, observations with high time 
resolution made with the Green Bank Telescope showed that the eclipse is mod¬ 
ulated at the spin period of pulsar Bp^ Modelling of the detailed eclipse profile 
by absorption in the doughnut-shaped closed-field-line region of the magnetosphere 
by Lyutikov and Thompson 84 allowed determination of the geometry of the binary 
system including the offset between the B-pulsar spin axis and the orbital angu¬ 
lar momentum axis, the so-called “misalignment angle” which they estimate to be 
about 60°. Even more remarkably, detailed measurements of the eclipse profile over 
about four years enabled Breton et al. 20 to directly estimate the rate of geodetic 
precession as 4?77 ± 0?66 yr _1 , consistent with the predicted precessional period of 
70.95 years based on GR (Equation 0). It is notable that no secular profile evolution 
has been observed for PSR J0737—3039A, implying that, unlike for the B pulsar, 
the misalignment angle for pulsar A is very small 

The most precisely measured post-Keplerian parameter is ui. This constraint is 
nearly orthogonal to that from the mass ratio R (Figure 0 giving values for the two 
neutron-star masses of mi = 1.3381 ± 0.0007 M 0 and m 2 = 1.2849 ± 0.0007 Mq.EI 
Together with the accurately measured Keplerian parameters, these two masses can 
be used to predict the remaining five post-Keplerian parameters using Equations 0 
0 As Figure0shows, GR gives a self-consistent description of the orbital motions, 
with all five independent constraints consistent with the masses derived from w and 
R. 

The most stringent test comes from the derived value of s, 0.99974 ± 0.00039. 
The ratio of the observed and predicted values of s is 0.99987 ± 0.0005, by far 
the most constraining test of GR in the strong-field regime. Furthermore, this test 
is qualitatively different to that for PSR B1913+16 in that it is based on a non- 
radiative prediction. With just 2.5 years of data analysed, the orbital decay term 
for the Double Pulsar is not measured as precisely as that for PSR B1913+16. 
However, since the phase offset grows quadratically (Figure 0 and the effect of 
receiver noise decreases as T 1 / 2 where T is the data span (assuming approximately 
uniform sampling), the precision of the Pf, measurement should improve as T 5 / 2 . 
Furthermore, J0737—3039A/B is much closer to the Sun than PSR B1913+16, thus 
reducing the magnitude and uncertainty of the correction due to differential Galactic 
acceleration. Also, using very long-baseline interferometry (VLBI), Deller et alP^ 
have shown that the transverse space velocity of the binary system is small, just 
24^g km s _1 , and so the Shklovskii correction to the observed Pf, is also small and 
accurately known. For these reasons, the orbit-decay test with the Double Pulsar 
will not be limited by these corrections, at least for another decade or so. Beyond 
that, improved models for the Galactic gravitational potential can be expected from 
analysis of GAIA data^^O anc j improved parallax and proper motion measurements 
may come from further VLBI observations, both reducing the uncertainty in these 
kinematic corrections. 

Several other post-Keplerian parameters can in principle be observed in binary 






April 21, 2015 0:29 WSPC/INSTRUCTION FILE grlOO'rnm'ijmpd 


14 



Mass of Pulsar A (M gun ) 

Fig. 8. Plot of companion mass (m 2 ) versus pulsar mass (mi) for PSR J0737—3039 with observed 
constraints interpreted in the framework of GR. The inset shows the central region at an expanded 
scale, illustrating that GR is consistent with all constraints. (M. Kramer, private communication) 


systems,^ testing other aspects of relativity. Kramer and We>P 5 show that only 
one of these, 5g, which quantifies relativistic deformations of the binary orbit, is 
potentially measurable with the Double Pulsar system. Even this will take more 
than a decade to reach an interesting level of precision given reasonable expectations 
about future observations. 

A potentially exciting prospect is to use measurement of higher-order terms for 
relativistic periastron precession to put constraints on the moment of inertia of a 
neutron star. 32 At the second post-Newtonian level, the periastron precession has 
three components: 

ill = h\pN + ll2pN + hso ( 9 ) 

where h \ v n is given by Equation 0 li 2 pN is the 2nd post-Newtonian (u/c) 4 con¬ 
tribution and hso gives the contribution from spin-orbit coupling. 16 32 This latter 
term depends on the moment of inertia and spin of the A pulsar (the spin angular 
momentum of the B pulsar is negligible) and can in principle be measured if the two 
leading components in Equation [9] can be measured to sufficient precision. Kramer 
and Wed 7 ® show that, given expected advances in precision timing and astrometry, 
a significant result could be obtained in 20 years or so. 

Spin-orbit coupling can also lead to a nonlinear variation in w as a function of 
time and a time variation in the projected semi-major axis of pulsar A, a\ sin 
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However, the nonlinear terms depend on sind, where 6 is the misalignment angle. 
As mentioned above, it seems as though 9 a is very small, so unfortunately these 
terms may be difficult to detect. 

2.1.4. Measured post-Keplerian Parameters 

Table [l] summarises the measured post-Keplerian parameters for the eleven pulsar 
binary systems where three or more post-Keplerian parameters have been mea¬ 
sured. Since only two measurements are required to determine the stellar masses, 
these systems provide the opportunity for tests of theories of relativistic gravity. 
Of the eleven systems, in six cases the companion star is believed to be a neu¬ 
tron star, in a further four cases a white dwarf companion is more likely and in 
one case the companion is probably a main sequence star. This latter system, PSR 
J1903+0327, evidently has had an unusual evolutionary history as a triple system in 
which one component was ejected.®® It is currently not clear if there is a significant 
contribution to the observed to from kinematic effects and/or spin-orbit coupling, so 
the utility of this system for tests of relativistic gravity is limited. There a further 
ten pulsar systems in the literature where just two post-Keplerian parameters have 
been measured, enabling estimates of the stellar masses, but no tests of relativistic 
gravity. 

Although an independent measurement of relativistic spin precession has been 
possible for just two systems so far, as shown in TableHJ secular changes in observed 
pulse profiles have been observed for four other pulsars and attributed to geodetic 
precession of the pulsar spin axis. 

2.2. Tests of Equivalence Principles and Alternative Theories of 
Gravitation 

Pulsars and especially binary pulsars have unique advantages in testing theories of 
relativistic gravitation as a result of their often rapid spin, short orbital periods and 
the ultra-high density of the under-lying neutron stars. As we have shown above, GR 
has been amazingly successful in describing all measurements to date. Never-the- 
less, investigations of quantum gravity and cosmology suggest that, in some regimes, 
extensions or modifications of GR may be required. This strongly motivates a search 
for departures from GR within existing experimental capabilities. 

Gravitational theories have Equivalence Principles at their heart. The Weak 
Equivalence Principle (WEP) is basic to Newtonian gravity, stating that acceler¬ 
ation in a gravitational field is independent of mass or composition. The Strong 
Equivalence Principle (SEP) adds Lorentz invariance (no preferred reference frame) 
and position invariance (no preferred location) for both gravitational and non- 
gravitational interactions. GR satisfies the SEP whereas other theories may violate 
the SEP or even the WEP in one or more respects. 

Comparison of different gravitational theories has been greatly facilitated by 
the “parametrized post-Newtonian” (PPN) formalism which describes observable 
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Table 1. Binary pulsars with three or more significant measured post-Keplerian parameters 


Pulsar/Parameter 

J0437—4715 

J0737—3039A/B 

J0751+1807 

J1141-6545 

B1534+12 

J1756—2251 

Peri, advance Co (°yr —1 ) 

0.016(8) 

16.8995(7) 

- 

5.3096(4) 

1.7557950(19) 

2.58240(4) 

Time dilation 7 (ms) 

- 

0.3856(26) 

- 

0.773(11) 

2.0708(5) 

0.001148(9) 

Orb.P deriv. P 5 (10 —12 ) 

3.73(6) b 

-1.252(17) 

-0.031(5) 

-0.403(25) 

-0.1366(5) 

-0.229(5) 

s = sin i 

0.6746(28) b 

0.99974(39) 

0.90(5) 

- 

0.9772(16) 

0.93(4) 

Comp, mass m 2 (Mq) 

0.254(18) 

1.2489(7) 

0.191(15) 

- 

1.35(5) 

1 . 6 ( 6 ) 

Geod. prec. £l p (°yr — - 1 ) 

- 

4.77(66) c 

- 

Note d 

0.59(10) 

- 

Binary companion a 

He WD 

NS 

He WD 

CO WD 

NS 

NS 

References 

□m 

Emu 

EMizl 

El 

[49] 

[46] 

Pulsar /Parameter 

J1807-2459B 

J1903+0327 

J1906+0746 

B1913+16 

B2127+11C 


Peri, advance Co (°yr —1 ) 

0.018339(4) 

0.0002400(2) b 

7.5841(5) 

4.226598(5) 

4.4644(1) 


Time dilation 7 (ms) 

26(14) 

- 

0.470(5) 

4.2992(8) 

4.78(4) 


Orb.P deriv. Pb (10 —12 ) 

- 

- 

-0.56(3) 

—2.423(l) b 

-3.96(5) 


s = sin i 

0.99715(20) 

0.9759(16) 

- 

- 

- 


Comp, mass m 2 (Mq) 

1.02(17) 

1.03(3) 

- 

- 

- 


Geod. prec. (°yr —1 ) 

- 

- 

Note d 

Note d 

Note d 


Binary companion 3. 

CO WD(?) 

MS 

NS(?) 

NS 

NS 


References 

E2 

[50] 

□miii] 

1145111431 

[62] 



Note : a: Binary companion types: CO WD: Carbon-Oxygen White Dwarf; He WD: Helium White Dwarf; MS: Main-sequence 
star; NS: Neutron star 

b: Dominated or significantly biased by kinematic effects 
c: For PSR J0737-3039B 

d: Effects of precession observed, but no independent determination of Q p 


or potentially observable phenomena in a theory-independent way. This formalism 
was first developed by Will and Nordtvedlj 150 for “weak-field” situations, that is, 
where e ~ GM/Rc 2 <C 1, where G is the Newtonian gravitational constant, M and 
R characterise the size and mass of the object and c is the velocity of light. Many 
gravitational tests are performed within the solar system where e < 10~ 5 , firmly 
in the weak-field regime. However in the vicinity of a neutron star e ~ 0.2 and so 
strong-field effects are potentially important. A number of authors have considered 
the generalisation of the PPN formalism to strong-field situations (see, e.g., Refs 
1311H4811351) allowing investigation of these effects. 

Some Lorentz-violating theories predict a dependence of the velocity of light on 
photon energy or polarisation. 149 Pulsar observations can be used to limit these the¬ 
ories, but potentially stronger limits come from observations of gamma-ray bursts, 
polarised extra-galactic radio sources and the cosmic microwave background. 95 

A recent comprehensive review of observational limits on theories of gravitation 
by C. M. Will can be found in Ref. 11491 Pulsar tests of relativistic gravitation have 
been reviewed by I. H. Staid 12 ® and, more recently, by N. WexP^H Further details 
on many of the topics discussed here may be found in these reviews. 
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2.2.1. Limits on PPN parameters 

The standard PPN formalism has ten parameters: 7ppn, /3ppn, £, a i, ct 2 , a 3 , ( 4 , 
£ 2 , C 3 and C-t- I n GR 7ppn, describing space curvature per unit mass, and /3ppn, 
describing superposition of gravitational fields, are unity and all others are zero. £ 
describes preferred location effects, the a n preferred frame effects and the others 
describe violations of momentum conservation. (03 also may be non-zero in this 
case.) Pulsar observations do not directly constrain 7ppn or /3ppn but are important 
in constraining many of the remaining PPN parameters and in fact currently place 
the strongest constraints on several parameters. 

Damour and Schafer 33 -! recognised that wide-orbit low-eccentricity binary pul¬ 
sars, which generally have a white-dwarf companion, could provide a valuable test 
of the SEP through a strong-field extension of the solar-system tests pioneered by 
NordtvedtP^ A violation of the SEP would cause bodies with different gravitational 
self-energy to fall at different rates in an external gravitational field. In a pulsar - 
white-dwarf binary system, this results in a forced eccentricity in the direction of 
the gravitational field, that of the Galaxy in this case. This eccentricity is given by 

_ 3 A §-L n n \ 

6F 2 wa(2ir/Pb) ( 0) 

where gj_ is the projection of the Galactic gravitational field on to the orbital plane, 
iL is the relativistic periastron advance, a is the orbit semi-major axis and Pb the 
orbital period. The ratio of the gravitational mass m g and the inertial mass m 
(which is exactly one if the SEP is obeyed) for body i is described by 

(—) = 1 + A, = 1 + r] N (E g /mc 2 )i + r]' N {E g /mc 2 ) 2 + ... (11) 

where ?pv is the Nordtvedt parameter, a function of several PPN parameters (see, 
e.g., Ref. 1149) and A = Ai — A 2 . Violations of the SEP will result in non-zero A. 

Since Damour and Schafer ® 3 first proposed this method of testing SEP violations, 
the number of suitable pulsar - white-dwarf binary systems has increased greatly. 
Gonzalez et alP^ combined data for 27 systems to place a 95% confidence upper 
limit on |A| of 4.6 x 10 -3 @ Since E g /mc 2 ~ 0.1 for a neutron star, this limit is 
not as strong as the weak-field limit on t]n < 3 x 10 -4 from lunar-laser ranging?® 
However, it does enter the strong-field regime and test possible violations of the 
SEP that solar-system tests cannot reach. 

The recent discovery by Ransom et alpAll 0 f a remarkable triple system contain¬ 
ing a pulsar, PSR J0337+1715, and two white-dwarf stars in essentially coplanar 
orbits, one in a relatively tight 1 . 6 -day orbit with the pulsar and the other in a 
much wider 327-day orbit around the inner system, opens up the possibility of a 
much more sensitive test of the SEP. Precise timing observations of the pulsar have 
already shown that the motion of the inner system is strongly affected by the grav¬ 
itational field of the outer white dwarf. The gravitational field of this star, which 


f This limit may be slightly under-estimated - see Wex. LAI 
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has a mass of about 0.41 M 0 , at the inner system is at least six orders of magni¬ 
tude larger than the Galactic gravitational field at the position of a typical pulsar 
and so the effect of SEP violations may be expected to be correspondingly largerf' 51 ^ 
Observations over several orbital periods of the outer star will almost certainly be 
necessary to isolate any SEP-related effects from other orbit perturbations. 

The wide-orbit low-eccentricity binary pulsars can also be used test for violations 
of local Lorentz invariance (LLI) of the gravitational interaction and momentum 
conservation that are described by the parameter a 3 and its strong-field generalisa¬ 
tion < 33 . Bell and Damoui 13 show that such violations produce a forced eccentricity 
analogous to that produced by the Nordtvedt effect given by 


|gf| = a3 


cpM p b ° 2 

24-7T P G(toi + m 2 ) 


sin/3 


( 12 ) 


where c p is a dimensionless “compactness” parameter, for neutron stars about 0 . 2 p^ 
and /3 is the angle between w, the velocity of the system with respect to a reference 
frame defined (for example) by the cosmic microwave background, and the pulsar’s 
spin axis. A similar analysis to that for the generalised Nordtvedt effect resulted in 
a 95% confidence limit of 1 0 : 3 1 < 4.0 x 10 -20 , some 13 orders of magnitude lower 
than the best solar-system testP^SIlt is worth noting that the observed small scatter 
in the period derivatives of MSPs had already been used by Bell^to limit |a 3 | to 
< 5 x 1CT 16 . 

Strong-field limits on the other two PPN parameters describing preferred frame 
effects, specifically LLI, di and & 2 , can also be obtained from low-eccentricity binary 
pulsar observationsNon-zero di induces a forced eccentricity in the direction of 
motion of the binary system velocity w, analogous to Nordtvedt d 3 tests discussed 
above, whereas non-zero 0:2 induces a precession of the orbital angular momen¬ 
tum about w. The best current limits come from observations of the binary pulsar 
systems PSRs J1012+5307 and J1738+0333, both of which have short orbital pe¬ 
riods, ~ 0.60 and ~ 0.35 days respectively, and extremely low orbital eccentricities 
e < 10 ~ 7 . 121 Furthermore, both of these pulsars have optically identified binary 
companions which, together with proper motion measurements from the timing 
observations, allow the three-dimensional space velocity of the binary system to 
be determined. For these pulsars, the observational data span is sufficiently long 
(iS 10 years) that relativistic periastron advance has significantly changed the ori¬ 
entation of the intrinsic eccentricity vector relative to the direction of any forced 
eccentricity, resulting in a potentially detectable change in the total eccentricity. 
The strongest limit on di comes from observations of PSR J1738+0333 as shown 
in Figure [HI conservatively a± < 4 x 10~ 5 . This is not only better than the best 
weak-field limit of a± < 2 x 10 -4 from lunar laser ranging®^ but also constrains 
strong-field effects as well. 

A precession of the orbital axis about the system velocity vector w induced 
by a non-zero value of 0:2 is potentially observable as a time variation x in the 
projected semi-major axis of the pulsar orbit x = a± sini. x is one of the possible 
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4 a 2 =0 i < 90° 



-4 


0 90 180 270 360 

Q. [ degree ] 


Fig. 9. Constraints on the strong-field PPN parameter &i from timing observations of the low- 
eccentricity binary pulsar PSR J1738+0333. The constraints are a function of the unknown ori¬ 
entation of the binary system on the sky (described by the longitude of the ascending node Q) 
and the unknown sign of sin i. If the PPN parameter 02 is assumed to be zero, then only certain 
ranges of Q are permitted. The shading corresponds to 68%, 95% and 99.7% confidence limits on 
di. CRef. ImTl 


post-Keplerian parameters in standard timing solutions and significant values have 
been measured for both PSRs J1012+5307 and J1738+0333! 121 There are several 
possible contributions to the observed x but in Ref. [1211 it is shown that all of these 
are negligible in these systems except that due to the changing orbit inclination i 
resulting from proper motion of the system. This is a function of the unknown Q 
values and for certain Q values there is no constraint. Consequently only a prob¬ 
abilistic limit on 0.2 can be obtained. By combining results for the two pulsars, a 
95% confidence limit of 1 0 : 2 1 < 1-8 x 10 -4 is obtained. 

This is not as constraining as a solar-system limit 1 0 : 2 1 < 2.4x 10~ 7 obtained from 
the present deviation of the Sun’s spin axis from the solar-systenr orbit normal," 
a limit that rests on the assumption that the two axes were aligned at the time of 
formation of the solar system. 

However, pulsars provide an even stronger constraint based on the stability of 
the spin axis of isolated pulsars. Any precession of the spin axis of a pulsar is likely to 
result in changes in the observed pulse profile (as observed with geodetic precession 
in binary pulsars as discussed in )J2] above). Sliao et al. 120 l compared mean pulse 
profiles for the isolated MSPs B1937+21 and J1744—1134 taken 10 - 12 years apart 
with the same observing system and found no perceptible change in the pulse width 
at 50% of the peak amplitude. They interpreted these results by assuming a circular 
beam profile for the main pulse in PSR B1937+21 and for PSR J1744—1134. All of 
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the angles in the problem can be estimated from modelling of radio and gamma- 
ray observations of the pulsar, taken together with known direction of the system 
velocity w with respect to the cosmic microwave background, except the angle 
of the projected spin axis on the sky. Probability histograms for ot 2 allowing for 
this unknown angle are shown in Figure flTIl The final result for the 95% confidence 
upper limit is I& 2 I < 1.6x 10 , about four orders of magnitude better than the limit 

described above based on orbital precession in pulsar binary systems and two orders 
of magnitude better than the limit based on solar spin precession. The assumption 
of circular beams in these MSPs is problematic, since there is good evidence for 
caustic enhancement in MSP pulse profiles, both radio and gamma-ray which 
would tend to elongate the beam in the latitude direction, reducing the effect of 
precession on the observed pulse profiles. However, taking this into account would 
probably increase the limit by a factor less than ten, and so this limit would remain 
the best available. 



-4 -2 0 2 4 

« 2 [ 10 ' 9 ] 


Fig. 10. Constraints on the strong-field PPN parameter a .2 based on the long-term stability of 
the mean pulse profiles for the isolated MSPs B1937+21 and J1744—1134. (Ref. [120) 

The stable pulse profiles of PSRs B1937+21 and J1744—1134 have also been 
used to limit the PPN parameter £ describing local position invariance (LPI), also 
known as the Whitehead parameter, and its strong-field counterpart £ (Ref. 1122D . 
The centripetal acceleration of Galactic rotation results in an anisotropy in the local 
gravitational field at a pulsar resulting in a precession of the pulsar spin vector 
around the direction of the Galactic acceleration with period. 



where vq is the velocity of Galactic rotation at the pulsar and ip is the (unknown) 
angle between the pulsar spin vector and the Galactic acceleration. Combining the 
results for the two pulsars in a way analogous to that for the 02 test described 
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above, Shao and Wexf23 obtain a 95% confidence limit ||| < 3.9 x 10 -9 . Even 
with the qualification about beam shapes mentioned above, this limit is at least 
two orders of magnitude better than the next best limit obtained from considering 
of the evolution of the solar spin-misalignment angle over the lifetime of the solar 
systemP^ 

The PPN parameter £2 is one of a number of parameters that may be non-zero 
in gravitational theories that violate conservation of total momentumA non¬ 
zero (2 would result in an acceleration of the centre of mass of a binary system 
that changes direction with periastron precession. This is best measured by looking 
for a change in the rate of orbital decay in a binary system with large periastron 
precession and a long data span. Wih 147 ^used observations of PSR B1913+16 with 
a 15-year data span to obtain a limit |£ 2 | <4x 10 -5 . Obviously, this test could be 
made more stringent using the long data spans now available for PSR B1913+16 
and other double-neutron-star systems with high ui. 


2.2.2. Dipolar gravitational waves and the constancy of G 

As described above, because of its tiny eccentricity, the pulsar - white-dwarf binary 
system PSR J1738+0333 has played a significant role in limiting the PPN parame¬ 
ters describing preferred frame effects. This system is composed of two very different 
stars, the neutron star we see as a pulsar and the companion which we know to be 
a white dwarf of mass 0.181 M 0 through its optical identification? 7 This and its 
short orbital period (~ 8.5 h) make it an ideal system for testing theories of gravity 
that predict a dipolar component to gravitational-wave damping as well as general 
scalar-tensor theories.^ 

From analysis of about 10 years of timing data obtained at Parkes and Arecibo 
observatories, Freire et al. 52 ^ find an observed rate of orbital decay for PSR 
J1738+0333 of P b = (—17.0 ± 3.1) x 10 -15 . To obtain the intrinsic rate of orbital 
decay, kinematic contributions from the differential acceleration of the binary sys¬ 
tem and the solar system in the Galactic gravitational field and the Shklovskii effect 
due to transverse motion must be subtracted, giving P b Int = (—25.9 ± 3.2) x 10~ 15 . 
Since the orbital parameters and the masses of the two stars are well known, the 
orbit decay due to GR can be accurately determined: P GR = (— 27.7f\'^) x 10~ 15 , 
leaving a residual orbit decay of P Res = ^.Olg'g) x 10 -15 . 

This residual orbit decay is consistent with zero, which can be interpreted as 
a further confirmation of the accuracy of GR. However, because of the very dif¬ 
ferent nature of the two stars in this binary system, this result also places strong 
constraints on theories of gravity that predict a dipolar component to gravitational- 
wave emission. Besides a dipolar component P b D , there are several other possible 
contributions to P b Int : 

pint = pM + pT + pD + pG (w) 

where P b M is due to mass loss from the binary system, P b T is a term resulting from 
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tidal effects on the white dwarf (tidal effects on the neutron star are negligible) and 
P b G is decay resulting from a possible variation in the gravitational ‘constant’ G.® 
Freire et al. 5 2l show that the likely M and tidal terms are small for this system, 
< 10 -15 , and so the limit on P b Int is effectively a limit on P b D + P G . 

Within certain restrictions on strong-field effects,® the dipole term is given by: 

= -^T Q m c ^- KD S 2 + G(s 3 ) (15) 

Pb <7+1 

where q = mi/rri 2 is the mass ratio (with subscript 1 refering to the neutron star), 
S = si — S 2 is the difference in “sensitivity” s of the mass of each body to a scalar 
field <j >, where 


/dlnm(</>) \ 

V dine/ ) ’ 


(16) 


and kd is a body-independent constant that describes the dipole self-gravity contri¬ 
bution in a given theory of gravity (see, e.g., Ref. 11491) . The sensitivity Si depends on 
the stellar equation of state and, for neutron stars, is typically about 0.15, whereas 
for a white dwarf it is ~ 10~ 4 . Therefore if kd is non-zero, dipole radiation will 
contribute to the orbit decay. 

The remaining term in the residual orbit decay is that due to possible variations 
in G. In weak gravity, G/G has been constrained to be less than 4 x 10~ 13 yr” 1 
from lunar laser ranging experiments,® giving 

P b G = -2^P b < 0.8 x 10" 15 . (17) 

Cr 

and hence \kd\ < 2 x 10 -4 . However, it is also possible to obtain independent esti¬ 
mates of the effects of dipole radiation and G by combining the J1738+0333 results 
with those for PSR J0437—4715.® This southern pulsar is in a wider orbit than PSR 
J1738+0333 and hence has a different mix of the dipole and G components, allow¬ 
ing them to be separated. After accounting for the fact that a changing G will also 
change the stellar massesp® the formal results are G/G = (—0.6+1.6) x 10 -12 yr -1 
and kjj = (—0.3 ±2.0) x 10~ 4 , both effectively upper limits. While the limit on kjj 
is the best available, the derived limit on G/G is about an order of magnitude 
weaker than the result (actually a limit on the variation of GMq) from the Mars 
Reconnaissance Orbited 1 and from lunar laser ranging.® 

Interestingly, pulsars have provided two other independent limits on G/G. 
ThorsettP® used determinations of neutron star masses from timing observations of 
double-neutron-star systems that formed many gigayears ago. In standard forma¬ 
tion scenarios, the mass of a neutron star depends on the Chandrasekhar mass, the 
maximum possible mass of a white dwarf star, just prior to the collapse to a neutron 
star. The Chandrasekhar mass is proportional to G~ 3 ^ 2 and so the observed small 
range of neutron star masses implies that G/G < 4 x 10 -12 . This limit has been 
somewhat weakened by recent discoveries of both less massive and more massive 
neutron stars in pulsar binary systems (see Ref. [70] for a recent review). 
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The very small observed rate of change of pulsar period P observed in some pul¬ 
sars (after correction for kinematic effects) may be used to set a further independent 
limit. 153 A variation in G will result in an inverse variation in the stellar moment 
of inertia, with the exact relation depending on the neutron-star structure. If the 
observed (intrinsic) P is entirely attributed to this effect, a limit of G/G < 2 x 10 -11 
is obtained. 

2.2.3. General scalar-tensor and scalar-vector-tensor theories 

Many alternate theories of gravity can be expressed in a “tensor-scalar” framework 
in which a scalar field <j> contributes to the “physical metric” through a coupling 
function A(<j>): 

g^ = A 2 {<t>) gill/ (18) 

where g is the usual tensor metric. The coupling constant may be expressed in 
different waysp*^ 0 ne of which is as an expansion around the asymptotic value of 
the scalar field cj) o: 

In A(<f>) = lnA((j> 0 ) + a o (0 - </>o) +/3 0 (^~ </>o) 2 + ••• (19) 

iRef. loTfl . For GR, ao = f3o = 0. In the well-known example of a tensor-scalar theory, 
the Brans-Dicke theory, the scalar coupling is described by a single parameter ujbd- 
For ujbd — > oo, the Brans-Dicke theory approaches GR. For this theory ao = 
1/(2ojbd + 3) and fio = 0. In other theories both the linear and quadratic terms in 
Equation [T9] (and higher-order terms) may be non-zero. 

The various observational constraints on PPN and post-Keplerian parameters 
can be expressed as limits in the (ao,fio) space for tensor-scalar theories as shown in 
Figure mi 5 - The Cassini experiment 14 placed a strong limit on the PPN parameter 
"/ppn = 1 + (2.1±2.3) x 10 -5 which translates to a limit on |ao| of about 0.003. Only 
the limits on dipole gravitational radiation from the asymmetric binary systems PSR 
J1141— 654d 4 ®and PSR J1738+033fP® rival the Cassini limit over most of the space. 
Since the precision of these measurements will increase with time, it seems likely 
that ultimately binary systems such as these will provide the strongest constraints 
on tensor-scalar theories. 

PSR J0348+0432 and its white-dwarf companion form another asymmetric bi¬ 
nary system, one that is distinguished by its very short orbital period (2.46 hours) 
and massive neutron star (2.01 ± 0.04 M 0 )P The high neutron-star mass is of in¬ 
terest since some scalar-tensor theories predict a strongly non-linear relationship 
between the strength of dipole GW emission and neutron-star mass or self-gravity 
(cf., Equation[l5|. For most of the parameter space of the class of theories discussed 
by Freire et al., 52 the limits on effective scalar coupling from PSR J0348+0432 are 
currently not as strong as those from PSR J1738+0333 but, because of the high 
neutron-star mass, they place stronger limits on some other theories with a greater 
degree of non-linear couplingP^H 
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Fig. 11. Constraints on the scalar-field parameters in the (|ao|,A)) plane from various observa¬ 
tional tests. Only the region below each line is allowed by the corresponding test. “SEP” refers 
to the test of the strong equivalence principle based on low-eccentricity pulsar - white-dwarf bi¬ 
nary systems“LLR” refers to lunar laser ranging resultand Cassini to the “Shapiro delay” 
experienced by signals to and from the Cassini spacecraft (on its way to Saturn) as its line of 
sight passed close to the Sun!^ Other binary-system tests are labelled according to the pulsar 
concerned. In GR, both ao and /3q are zero. (Ref. 152b 


BekensteirPS has proposed a relativistic generalisation of the so-called “MOND” 
theory of gravity 91 that seeks to avoid the need for dark matter in galactic dynamics. 
The generalisation invokes an additional vector field and hence is known as a tensor- 
vector-scalar theory. Such theories relax some of the constraints on tensor-scalar 
theories, in particular, the dipole radiation constraints, and allow significantly larger 
values of ao- However, as shown by Freire et al.pl the binary pulsar results still 
significantly constrain theories of this type and in fact are more constraining than 
solar-system tests. With future observations, binary pulsar tests have the potential 
to make this class of theories untenable. 

In another example of the use of pulsar observations, especially the limits on 
dipolar-GW radiation, to place limits on gravitational theories, Yagi et al. 155 have 
strongly limited the allowed parameter space for the LLI-violating “Einstein-iEther” 
and the “Khronometric” theories. 
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2.3. Future Prospects 

As described above, GR has provided an accurate description of all pulsar timing 
results obtained so far. However, continued refinement of existing methods and de¬ 
velopment of new tests is highly desirable. Continued pulsar timing measurements, 
especially with the advent of new and more sensitive observing facilities such as 
the 500-rn Arecibo-type FAST radio telescope in China^and the Square Kilometre 
Array (SKA) in South Africa and Australia^ will certainly improve on existing 
limits and enable new tests of gravitational theories. They may even demonstrate 
a failure of GR and hence a need for a modified or conceptually different theory of 
gravity. Conversely, if GR is assumed to be valid, results of astrophysical significance 
can be deduced from the observations. For example, as described in Section 12.1.31 
observations of higher-order terms in the relativistic perturbations may enable a 
measurement of the moment of inertia of a neutron star. 

These more sensitive radio telescopes can also be used to search for previously 
unknown pulsars and binary systems that are suitable for tests of gravitational 
theories. Past experience has shown that pulsar searches repeatedly turn up new 
classes of object. This potential is wonderfully illustrated by the recent discovery of 
the pulsar triple system PSR J0337+1715 which promises to provide a strong limit 
on violations of the Strong Equivalence Principle. A dream for such searches is the 
discovery of a pulsar in a close orbit around a black hole as this would offer much 
more stringent tests of gravity in the strong-field regime. 80 Discovery of a pulsar 
orbiting the black hole at the centre of our Galaxy with an orbital period of a few 
months or less could even allow a test of the so-called “no-hair” theorem for black 
holes PD 

3. The Quest for Gravitational-Wave Detection 

The direct detection of the gravitational waves (GWs) predicted by Einstein’s gen¬ 
eral theory of relativity and other relativistic theories of gravity is one of the major 
goals of current astrophysics. As described in $2] above, we have excellent evidence 
from the orbital decay of binary systems for the existence of GWs at the level pre¬ 
dicted by GR, but up to now there has been no direct detection of the changing 
curvature of spacetime induced by a passing GW. This changing curvature induces 
a change in the proper distance between two test masses, described by the gravi¬ 
tational strain h = 8L/L. The problem is that, for any likely source, h is tiny. For 
example, the LIGO gravitational-wave detector 2 hopes to detect the merger of two 
neutron stars at a distance of 100 Mpc for which h ~ 10 —22 , a change in the length 
of its 4-km arms of 10~ 18 nr or 10~ 3 of the diameter of a proton! 

Pulsar timing can measure a change in the proper distance between the pul¬ 
sar and the telescope. Systematic changes in timing residuals for a given pulsar 
reflect unmodelled changes in the effective time of emission, the pulsar position, 
the propagation path or the position of the telescope. With care, and with observa¬ 
tions of multiple pulsars, residual delays due to changing proper distances can be 
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isolated, effectively giving us a set of interferometers with baselines of > 10 16 km! 
However, even in the best cases, we can only measure the interferometer “phase” 
to about 100 ns, so that the limiting strain is about 10 -18 . Unlike ground-based 
laser-interferometer systems, which are most sensitive to GW signals with frequen¬ 
cies around 100 Hz, pulsar timing systems are most sensitive to signals with fre¬ 
quencies around the inverse of the data span, typically a few nanohertz. Potential 
sources of detectable GWs in this low-frequency band include super-massive black- 
hole (SMBH) binary systems in distant galaxies and cosmic strings in the early 
Universe. 


3 . 1 . Pulsar Timing Arrays 

The effect of GWs on pulsar timing signals was first considered by SazhirP^ and 
DetweilerJ^ with the latter being the first to consider the effect of a GW from a 
distant source passing over a pulsar and the Earth. For this case, it can be shown 
that the net effect on the observed pulsar arrival times is simply the difference 
between the effect of the GW passing over the pulsar and the effect of the GW 
passing over the Earth. For a GW travelling in the z direction, the redshift z of the 
pulse frequency v for a pulsar at distance d with direction cosines (a, /3, 7) is given 
by: 


z{t) 


vp ~ v{t) 

Vp 


a 2 — /3 2 
2(1 + 7) 


Ah + + 


a/3 
1 + 7 


Ah x 


( 20 ) 


where Ah a = h p A — h A , with A representing the two possible wave polarisation 
states (+,x) in GR, and h p A (t — djc ) and h A (t) the gravitational strain at the 
pulsar and the Earth, respectively]^ The observed timing residuals are then given 
by the integral of the redshift, 

R(t) = f z(t')dt'. (21) 

J 0 

The coefficients multiplying the A h + and A h x terms are the “antenna patterns” 
for the two polarisations as illustrated in Figure [12] A pulsar located in the GW 
propagation direction has zero sensitivity to the GW since its direction cosines a and 
/3 are zero. Furthermore, despite the (1+7) term in the denominator of Equation [20] 
the response for a pulsar exactly in the —z direction (the same direction as the GW 
source) is also zero. This comes about because of a cancellation of the (1 + 7) term 
by the expansion of Ah a when 1 + 7 <C 1 (cf. Ref. [78]). 

A pulsar timing array (PTA) consists of a set of pulsars spread across the sky 
which have precise timing measurements over a long data span. The detection of 
GWs by PTAs depends on the correlated timing residuals for different pulsars given 


8 See Ref. [5] for a rederivation of the Detweiler result. 
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Fig. 12. Effective “antenna pattern” for detection of a GW with pulsar timing. The wave is 
propagating in the +z direction and is assumed to have the + polarisation. The pattern for the x 
polarisation is the same but rotated by 45° about the 2 : axis. (Ref. I25D 


by the Earth term h^(t) in Equation[20] GWs passing over the pulsars produce un¬ 
correlated residuals because of both the retarded time and the different GW envi¬ 
ronment for each pulsar. Also the pulsars themselves have uncorrelated timing noise 
at some level, either intrinsic or resulting from uncorrected variations in interstel¬ 
lar delays. Because the expected GW strain is so weak, only MSPs have sufficient 
timing precision to make GW detection with PTAs feasible. 

For an isolated source of continuous GWs, say an SMBH binary system in a 
nearby galaxy, in principle, both the Earth term and and the pulsar term in Equa¬ 
tion [20] could be detected. For a rapidly evolving source, the pulsar term and the 
Earth term may be at different frequencies because of the retarded time of the pul¬ 
sar term (see, e.g., Ref. [65l)- They can then be added incoherently to increase the 
detection sensitivity. For a non-evolving source, i.e., 5fb <C 1/T where Sfb is the 
change in binary orbital frequency over the pulsar timing data span T, in principle 
the pulsar term and the Earth term could be summed coherently for optimal sensi¬ 
tivity. As is discussed further in Section [3.41 below, unfortunately we currently don’t 
know enough pulsar distances to sufficient accuracy to make this coherent addition 
possible. In a PTA, the pulsar terms therefore add with random phase, washing out 
the fringes in the antenna pattern (see also Ref. 1781) and adding “self-noise” to the 
signal from the Earth term. Since the antenna pattern (Figure [T2]) has a maximum 
for pulsars roughly in the same direction as the GW source, the maximum response 
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of a PTA is toward the greatest concentration of pulsars in the array. 

A stochastic background of nanohertz GW from many SMBH binary systems 
in distant galaxies is likely to be the signal first detected by PTAs. To a first 
approximation, this background is also likely to be statistically isotropic, i.e., the 
expectation value ( h 2 ) is independent of direction when averaged over typical data 
spans. Hellings and DownJ^were the first to show that, in this case, the correlation 
between GW-induced timing residuals for two pulsars separated by an angle 9 on 
the sky is dependent only on 9 and not on the sky positions of the two pulsars. The 
zero-lag correlation function, commonly known as the Hellings & Downs curve and 
obtained by integrating the product of the antenna patterns (Figure Il2l) for the two 
pulsars over all possible GW propagation directions, is given by: 

1 3x (, 1 \ 

Chd = 2 " 2 ~ \ 11 X ~ 6 ) ’ (22) 

where x = (1 — cos#)/2, and is plotted in Figure fl3l chd goes negative for angular 
separations around 90° and then positive again for pulsars that are more-or-less 
opposite on the sky - this is a direct consequence of the quadrupolar nature of 
GWs. It is also important to note that the limiting value as 9 —> 0 is 0.5, not 
1.0. This is a consequence of the fact that the pulsar terms in Equation [JH] are 
uncorrelated and, on average, of equal amplitude to the Earth term. The scatter 
in the simulated correlations results from the random phases of the pulsar terms 
and illustrates the “self-noise” that limits the sensitivity of PTA experiments in the 
strong-signal limit. 


3 . 2 . Nanohertz Gravitational-Wave Sources 
3.2.1. Massive black-hole binary systems 

There is good evidence that massive black holes form in the centre of galaxies at 
very early times (see, e.g., Ref. 1140 1 and also that merger events play a major role 
in galaxy growth (see, e.g., Ref. 1132 ). When two galaxies, each containing a cen¬ 
tral massive black hole, merge, dynamical friction will result in the two black holes 
migrating to the centre of the merged galaxy to form a binary system, with an 
estimated timescale for the migration of the order of giga-years (see, e.g., Ref. l68l) . 
When the binary separation is less than about 1 pc, loss of energy to GWs becomes 
the dominant orbital decay mechanism and the binary system will ultimately coa¬ 
lesce, becoming a strong GW source as it spirals in. There remains much controversy 
about the efficiency of orbital decay mechanisms as the binary separation approaches 
a parsec - known as the “last parsec problem”. Some (see, e.g., Ref. l69 | Il07 l) ar¬ 
gue that dissipation mechanisms will quickly move the binary system through this 
phase, whereas others (see, e.g., Ref. (271 argue that the binary is likely to stall at 
separations where the gravitational decay is ineffective. Detection of a stochastic 
GW background (GWB) would resolve this issue. 

Since large numbers of binary systems with different orbital periods contribute to 
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Angular separation (deg) 

Fig. 13. The Hellings & Downs correlation function, i.e., the correlation between timing residuals 
for pairs of pulsars as a function of their angular separation for an isotropic stochastic background 
of GWs. Also shown are simulated correlations between the 20 pulsars of the Parkes PTA for a 
single realisation of a strong GW signal that dominates all other noise contributions. (Ref. I58D 

the GWB, it is a broadband signal which is best described in the spectral domain. It 
is convenient to express the amplitude of the GW signal in terms of the dimensionless 
“characteristic strain”, defined by 

h c = 2f\h(f)\ (23) 

where h(f) is the Fourier transform of h(t) and / is the GW frequency. (Note that, 
for a circular binary system, the frequency of the emitted GW is twice the binary 
orbital frequency / = 2/{,.) Two other quantities that are often used to parameterise 
GW spectra and detector spectral sensitivities are the square root of the one-sided 
strain power spectral density 

Sl /2 (f ) = hj- 112 (24) 

and the GW energy density as a fraction of the closure energy density of the Universe 

2tt 2 

Qgw = jTrfhlif) (25) 

where Hq is the Hubble constant. 

In order to understand the astrophysical implications of results obtained from 
PTA experiments, it is necessary to have estimates of the likely strength of signals 
from potential sources of nanohertz GW. For a cosmological population of SMBH 
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binary systems at luminosity distance Dl and redshift z, the local energy density 
in GW at frequency / is given by: 


POO PC 

fS E (f) = dz 

Jo Jo 


dMr 


drn 


1 


1 dE a 


dzdM c (1 + z) D 2 l dln/ r 


(26) 


where M c = (M 1 M 2) 3 / 5 (Mi + M 2) -1 / 5 is the binary chirp mass and Mi and M 2 
are the masses of the binary components, d 2 n/(dzdM c ) is the comoving density 
of binary systems with redshift and chirp mass between z and z + dz and M c and 
M c +dM c , respectively, and dE g /d In f r is the total energy emitted by a single binary 
system in the logarithmic frequency interval dlnf r , where f r = /(I + z) EMEHM] 
The local GW energy density is related to the local characteristic strain by 

(27) 


fS E (f) = 


and for a circular binary system 


d In f r 


G'2/ 3 7r 2 /3 


M 3 / 3 / 2 / 3 . 


Therefore we have 


Mr, = 

CU> 37tV3 C ! / 


dz 


dMr 


1 


1 


(28) 


■M c 5/3 . (29) 


'dzdM c (1 + z) 1 / 3 D\ 

As PliinnejhAd has emphasised, the result that h c oc f~ 2 / 3 for a cosmological 
population of circular binary systems decaying through GW emission is quite general 
and independent of any particular cosmology, black-hole mass function or galaxy 
merger scenario. Consequently the spectrum of the GWB is often parameterised as 
follows: 


hc(f) = A lyr (J-^j 


(30) 


where /i yr = (lyr)^ 1 , Ai yr is the characteristic strain at /i yr and a = —2/3 for the 
case described above. For pulsar timing experiments, the one-sided power spectrum 
of the timing residuals is given by 

p U) = ( 31 ) 

Consequently, a GWB produces a very “red” modulation of the timing residuals 
with a spectral index of —13/3 for a = —2/3. 

In order to estimate the likely strength of this modulation, the factor 
d 2 n/(dzdM c ) in Equation 1291 must be evaluated. This requires a prescription for 
the cosmological evolution of massive black-hole binary systems in galaxies. Dif¬ 
ferent approaches to this problem have been taken by different authors. An early 
paper by Jaffe and Backei®2l used observational constraints on close galaxy pairs 
coupled with a black-hole mass function, whereas another early paper by Wyithe 
and LoetP®] used a prescription for merger of dark-matter halos coupled with differ¬ 
ent scenarios for growth of massive black holes in galaxies. The latter approach was 
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developed further by Sesana et al. 115 who showed that the GWB spectrum steepens 
at frequencies above about 10 -8 Hz since the number of binary systems contribut¬ 
ing to the background at these frequencies becomes small. This is illustrated in the 
left panel of Figure [TT] which shows that binary systems at z < 2 contribute most 
of the strain to the GWB. The right panel shows that massive binary systems with 
M c > 10 8 M 0 contribute most of the low-frequency GWB but that only systems 
with M c ^ 10 8 Mg contribute to the high-frequency end. Importantly, at the high- 
mass end, only a few systems contribute to the GWB. As shown in Figure [HI this 
leads to a steepening and large uncertainties in the expected GWB spectrum at 
frequencies > ICG 8 Hz. These results were confirmed and extended by Sesana et 
aipUl arid Ravi et al using the Millenium simulation of the cosmological evolu¬ 
tion of dark matter structures 1281 to define a merger history together with various 
prescriptions for galaxy and black-hole formation and growth. 



redshift 



10 7 10 8 10 ® 
chirp mass [M 0 ] 


Fig. 14. Left panel: Number of binary black-hole systems and their contribution to the charac¬ 
teristic strain h c of the GWB as function of source redshift z. The solid lines are results from 
a Monte Carlo approach and the dotted lines are from a semi-analytic analysis. In both panels, 
the upper histograms are for a GW frequency / = 8 X 10“ 9 Hz and the lower histograms for 
/ = 10~‘ Hz. Right panel: The lower two panels are as for the left panel but as a function of 
source chirp mass M c . The upper panel shows the number of frequency bins spanned by the chirp 
over the 5-year span of the simulation. (Ref. 11151) 


While there is some consensus about the form of the GWB spectrum, there 
remain significant uncertainties. For example, Ravi et al. 10 ^ consider the effects of 
the stellar environment on the late evolution of massive black-hole binary systems in 
the cores of galaxies and conclude that the effect of dynamical friction is important, 
both in extracting energy from the binary system and inducing an eccentricity. Both 
of these have the effect of reducing the strength of the predicted GWB, especially at 
the low-frequency end, consequently making its detection by PTAs more difficult. 
On the other hand, McWilliams et alPS argue for a model in which all black- 
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hole growth is by merger rather than by accretion after coalescence, which is the 
main contributor to black-hole growth in the models discussed above. This leads to 
predictions of a significantly larger GWB characteristic strain compared to previous 
predictions and, hence, the imminent detection of the GWB by PTAs. 



observed frequency [Hz] 


Fig. 15. Characteristic strain spectrum for the GWB based on different prescriptions for black- 
hole growth by accretion between mergers (solid, dashed and dotted lines) and different black-hole 
mass functions (different panels). (Ref. fU6ll 


As Figure [15] indicates, there is a possibility that the nearby universe could con¬ 
tain a massive black-hole binary system with an orbital period of the order of a few 
years that actually dominates the nanohertz GW spectrum. This opens up the excit¬ 
ing prospect of the GW detection and study of an isolated supermassive black-hole 
binary system using pulsar timing and even the possibility of detection and study of 
the system in the electromagnetic bands - so-called “multi-messenger” astronomy. 
Searches for binary GW sources will be described in 1 )3.31 and the prospects for their 
detailed study will be discussed in 1 )3.41 

3.2.2. Cosmic strings and the early Universe 

Cosmic strings and the related cosmic super-strings are one-dinrensional topological 
defects which may have formed in phase transitions in the early Universe. Cosmic 
strings occur in standard field-theory inflation models, whereas superstrings are 
found in brane inflationary models. The idea that such strings will oscillate and 
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hence emit GWs was first proposed by Vilenkin Such oscillations may contribute 
to the stochastic GWB (see, e.g., Ref. 1251) or generate bursts of GW radiation from 
string cusps and kinks (see, e.g., Ref. |36j). The amplitude of GWs from cosmic 
strings is dependent on a large number of poorly known (or unknown) parameters 
and hence is very uncertain (see, e.g., Ref. I lllil . Key parameters are the string 
tension /x, usually parameterised by the dimensionless quantity G/x/c 2 , and the size 
a of string loops relative to the horizon radius at the time of birth. Other significant 
parameters for the GWB are the intrinsic spectral index q of the GW emission, a 
characteristic node number n* for the high-frequency cutoff in the emission spectrum 
and the probability p of “intercommutation”, that is, intersecting strings dividing 
and the two parts exchanging. Such intercommutation can, for example, form two 
smaller loops from an intersecting twist in a larger loop. For standard strings, p = 1 
but it may be less for superstrings. 

Vibrating cosmic strings are likely to decay by emission of GWs in a series of 
harmonics with fundamental frequency 2c/1, where l is the length of the loop, with 
an initial value aDn, where Dh is the horizon distance at the time of loop creation. 
The rate of energy loss for vibration mode n of a given loop is: 


dEgw 

dt 


7.-9 


= r 


e:=i «- 9 


G{i 2 c 


(32) 


where T is factor depending on the shape of the loop, typically about 50. 111 1 As the 
loop loses energy, it shrinks and eventually disappears. The creation of loops through 
intercommutation and their decay through GW emission sets up an equilibrium 
distribution of loop sizes. Sanidas et al. 111 compute the number density of loops 
as a function of loop length and time and hence, using Equation [32j the predicted 
spectrum of the GWB from string loops as a function of the various parameters. As 
Figure flfil shows, the spectrum is very broad, extending all the way from nanohertz 
to Megahertz. Since the lowest (n = 1) frequency is proportional to the string length, 
the weaker and smaller loops do not contribute to the nanohertz background. 


3.2.3. Transient or Burst GW sources 

The prime target of ground-based laser-interferometer GW detectors is the burst 
emitted at the coalescence of a double-neutron-star system. Such a burst is intense 
for just a few milliseconds and clearly cannot be detected by PTA experiments 
which typically are sensitive to signals of duration between a few weeks and a few 
years. Possible sources of longer-duration bursts are coalescence of SMBH binary 
systems, highly eccentric massive black-hole binary systems, and formation and 
decay of cusps and kinks in cosmic strings. 

A major difference between detection of GW bursts and continuous GW sources 
is that there is generally no interference between the Earth term and the pulsar 
terms in the signal detected by PTAs (Equation [20li . This is because the duration 
of the burst (by definition, less than the observational data span) is much less than 
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log io f (Hz) 

Fig. 16. Energy density spectrum for GWs emitted by cosmic strings as a function of the di¬ 
mensionless string tension Gfi/c 2 . Other parameters are held fixed at values a = 10 7 . q = 4/3, 
n* = 1 and p = 1. The dashed spectrum is for the critical point where TGfijc 2 at a; spectra above 
this are for large loops and spectra below are for small loops. (Ref. 11111) 


the light-time to the pulsars and so, except for a source in the same direction as a 
pulsar, the burst will occur at very different times in the Earth and pulsar terms. 
The pulsar term reflects the effect of the burst on the pulsar at a time d/c before the 
burst arrives at the Earth, but the pulsar term is detected at a time (d/c)(l + cos$) 
later than the Earth term, where d is the pulsar distance and 9 is the angle between 
the pulsar and the GW propagation direction as seen from the Earth. 

Figure [17] shows the GW waveform produced by the coalescence of two black 
holes in a coordinate system where all the signal is in h + . The maximum amplitude 
of the waveform is about 0.1 in the time units of Figure [lT] or about 0.1 cMTq/D 
or ~ 10~ 14 Mg/D g P c where Mg = 10 ~ 9 M is the total system mass in solar units 
and Dg P c is the (comoving) distance in Gpc. Although this is comparable to the 
strain sensitivities achieved by current PTAs for continuous GW signals (see T1.3I 
below), even for the largest SMBHs, the timescale of the burst, ~ 200 MTq is only 
of order 10 days and the period of the oscillation is about an order of magnitude 
less than that. Not only is this too short to be resolved by any existing PTA, but 
the sensitivity over this short interval would be much less than that achieved for 
CW signals integrated over the entire data span. Space-based laser interferometer 
systems such as the planned eLISA 3 have the potential to directly detect these 
bursts. 

However, Figure [T7] also shows another effect known as GW “memory” which 
is potentially detectable with pulsar timing. During the coalescence event, a non- 
oscillatory component to the gravitational strain builds up, so that at the end of 
the “ring-down” phase, the strain has a permanent offset from the pre-coalescence 
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Fig. 17. Gravitational waveform resulting from the coalescence of two equal-mass black holes 
(reduced mass ratio rj = 0.25) with total mass M at distance R. Both M and R are expressed 
in time units; the conversions to conventional units are M —> MTq and R —> R/c. 0 and are 
assumed source directions. The dashed line is the predicted waveform if the gravitational memory 
effect is ignored. (Ref. [44} 


value. The amplitude of the memory effect is 


_ ctjMTq 
~ D/c 


10- 15 


Mg 

D 9 


(33) 


where e ss 0.07 is the mass fraction contributing to the memory effect and 77 is the 
reduced mass fraction, 0.25 for equal-mass binary components UAlso] 

This step change in h produces a step change in the observed pulse frequency 
A v/v = h m , be., a “glitch”. This glitch persists until it is reversed by the pulsar 
term at a time (d/c)(l +cos 9) later. In a PTA, these reversals will occur at different 
times for different pulsars. Of course, it is also possible that a GW memory jump 
could be detected in the pulsar terms, but there it may be confused with a real glitch 
in the intrinsic pulse frequency, whereas in the Earth term there is a correlation in 
the effect on different pulsars. However, glitches in MSPs are rare (only one very 
small glitch detected so far: Ref. [55]). Also, real pulsar glitches are generally spin- 
ups, whereas a GW-memory jump may be of either sign, depending on pulsar-source 
angle. Therefore, as Cordes and .JenetP^ have discussed, the pulsar terms may give 
an improved probability of detection. 

Black-hole binary systems with circular orbits emit GW at the second harmonic 
of the orbital frequency, i.e., / = 2/&. For eccentric orbits, the GW emission be¬ 
comes more burst-like as the accelerations and hence GW power are greatest around 
periastron when the two black holes are closest together!^ i n the spectral domain, 
power spreads to higher harmonics and also to the fundamental frequency fb . In the 
gravitational-decay phase of evolution, when energy loss is dominated by GW emis- 
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sion, the orbit tends to circulariseP However, at earlier phases of the orbital decay 
when three-body stellar interactions or interaction with a gaseous disk surround¬ 
ing the binary system are important, the eccentricity may grow. Stellar three-body 
interactions can result in orbital decay through dynamical friction but probably 
result in a modest increase in the eccentricity of the black-hole binary system (e.g., 
Ref. [501) . However, Roedig et alP^find that initially mildly eccentric binary systems 
decaying through interaction with a gaseous disk evolve toward a limiting eccentric¬ 
ity in the 0.6 0.8 range. In a regime of frequent mergers it is even possible that a 

black-hole triplet could forrrP and, in this case, eccentricities as high as 0.99 could 
exist. 

Finn and Lomrnen 47 have investigated the GW emission and the resulting timing 
residuals from a close parabolic encounter of two massive black holes. As Figure UHl 
shows, an encounter of two 1O 9 M 0 black holes located at a distance of 15 Mpc 
with a minimum separation of 0.02 pc produces a GW burst of duration about 1 
year and maximum GW strain ~ 10~ 13 . This results in potentially detectable PTA 
timing residuals of about 1/zs amplitude with the same timescale. Unfortunately, 
the probability of having such a close encounter of two SMBHs in the local universe 
within PTA observational data spans is not high. 




Fig. 18. The upper plot shows the gravitational waveforms in the h -(_ and h- polarisations re¬ 
sulting from a parabolic encounter with an impact parameter of 0.02 pc of two 1 O 9 M 0 black holes 
located at a distance of 15 Mpc in the direction of the Virgo Cluster. The lower plot shows the 
resulting timing residuals for several PTA pulsars. (Ref. S3 


Cosmic strings are another potential source of GW burst emission. They can 
radiate over a wide frequency range with a huge range of possible amplitudes 
and timescales (cf., Figure fl6l) depending on the detailed mechanism invoked 
(loops, cusps, short strings, etc.) and the very wide (virtually unlimited) param¬ 
eter space.®- 27 77 Short bursts radiating in the LIGO and eLISA bands may be 
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frequent and unresolved, producing a GWB at these frequencies. However, bursts 
with longer timescales, producing radiation in the nanohertz band, are also possible 
but are likely to be extremely rareP^ Consequently, while in principle such bursts 
could be detected, in practice it is unlikely that PTAs will be able to significantly 
constrain models for GW burst emission from cosmic strings. 

3 . 3 . Pulsar Timing Arrays and Current Results 

In this section we first describe the three main PTAs currently operating world-wide: 
the European Pulsar Timing Array ( EPTA ), the North American pulsar timing 
array ( NANOGrav) and the Parkes Pulsar Timing Array ( PPTA ), and the collab¬ 
oration between them, the International Pulsar Timing Array ( IPTA ). PTAs have 
many possible applications such as establishing a pulsar-based timescale, 3 ' investi¬ 
gating the accuracy of solar-system ephenreridesp 6 and investigating the properties 
of the pulsars themselves (e.g., Refs [156 1 Ill8ft and of the intervening interstellar 
medium (e.g., Ref. E3- However, here we concentrate on what is undoubtedly their 
primary scientific goal, the direct detection of gravitational waves. Unfortunately, 
in common with other GW detection efforts around the world, PTAs have so far 
only been able to place limits on the strength of signals from potential GW sources. 
However, these limits are now beginning to seriously constrain the astropliysical 
source models and the assumptions that go into them and hence have implications 
that go far beyond the GW studies themselves. 

3.3.1. Existing PTAs 

The EPTA uses five large radio-telescopes in Europe, the Effelsberg 100-rn tele¬ 
scope in Germany, the Nangay Radio Telescope in France (95-rn equivalent area), 
the Westerbork Synthesis Radio Telescope in the Netherlands (similar effective area 
to the Nangay telescope), the 76-m Lovell Telescope at Jodrell Bank in England and 
the recently completed 64-m Sardinia Radio Telescope in Italy, to observe about 40 
MSPs with a cadence of between a few days and 30 days for different pulsars?- 3 Dif¬ 
ferent telescopes observe at different frequencies in the range 0.3 - 2.6 GHz, but all 
are instrumented at 1.4 GHz. Normally the five telescopes observe independently, 
but in a project known as the “Large European Array for Pulsars” (LEAP) 1.4 
GHz signals over a bandwidth of 128 MHz from the five telescopes can be summed 
coherently to form a 194-m equivalent diameter radio telescope. The different tele¬ 
scopes use different signal-processing systems, either digital filterbanks or coherent 
dedisperion systems or both. 

NANOGrav makes use of the 300-m Arecibo radio telescope in Puerto Rico and 
the 100-m Green Bank Telescope ( GBT) in West Virginia. 8 ' A sample of about 36 
pulsars is observed, typically at 3-week intervals. At Arecibo, observations are made 
in bands centred at 430 MHz and 1410 MHz, whereas at the GBT , the observed 
bands are centred at 820 MHz and 1500 MHz. Currently both radio telescopes use 
coherent dedispersion systems with bandwidths up to 800 MHz, but in the past a 
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range of filterbank and coherent dedispersion systems with more limited bandwidths 
have been used. 

As the name suggests, the PPTA uses the Parkes 64-m radio telescope located in 
New South Wales, Australia. A sample of 22 pulsars is currently being observed with 
regular observations at 2 - 3 week intervals in three bands around 730 MHz, 1400 
MHz and 3100 MHz respectively! 86 ^ 56 ^ Coherent dedispersion systems are used at 
730 MHz and 1400 MHz with bandwidths up to 310 MHz and digital filterbanks at 
1400 MHz and 3100 MHz with bandwidths of 256 MHz and 1024 MHz respectively. 

Data sets from all three PTAs have spans ranging from a few years up to about 
20 years for different pulsars; three (including the original MSP, PSR B1937+21) 
have Arecibo data spans of nearly 30 years. The three PTAs together observe about 
50 pulsars with some being observed by two or even all three of the PTAs. Their 
distribution on the sky is shown in Figure flT)l 



RA= 12 h 


Fig. 19. Distribution on the sky of MSPs being timed by the three PTAs, with different symbols 
for each PTA. Right ascension increases to the left with 0 h at the plot centre. The symbol size is 
related to the ratio S 1400 /-P> where S 1400 and P are the pulsar 1400 MHz flux density and pulse 
period respectively. (Ref. [85} 

Given that the combined data set of the three PTAs contains a larger number of 
pulsars, improved observation cadence and greater frequency diversity than the data 
set of any one PTA, there is a strong motivation to combine all the available data 
sets to obtain maximum sensitivity for PTA scientific objectives. The International 
Pulsar Timing Array ( IPTA ) consortium was set up to facilitate progress toward 
this goalP^The IPTA also arranges annual science meetings and student workshops 
and provides a forum for outreach programs and other activities related to PTA 
research. 

3.3.2. Limits on the nanohertz GW background 

As discussed in Section 13.2.11 GWs from a cosmological distribution of SMBH bi¬ 
nary systems are expected to contribute a very “red” signal to the spectrum of 
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pulsar timing residuals. The expected signal from other GWB sources is similar. 
Consequently, long-term observations of a single pulsar with little or no detectable 
intrinsic timing irregularities can be used to place a limit on the strength of the 
GWB in the Galaxy. Of course, statistical limits can be improved by using data 
from several such pulsars. An early limit on the GWB at a frequency of 4.5 nHz 
was set by Kaspi et al.® using Arecibo observations of two MSPs, PSR B1855+09 
and PSR B1937+21, with a 95% confidence limit on ttcwh 2 of 6 x 10 -8 , where 
h = Hq/100 km s -1 . 

Since the advent of the various PTA projects, both the quality and quantity 
of timing data sets has improved and a variety of analysis techniques have been 
employed to extract increasingly restrictive limits. Based on early PPTA data on 
seven pulsars combined with the Kaspi et al. PSR B1855+09 data set, Jenet et 
al.® used a “frequentist” approach with a statistic based on the amplitude of the 
low-frequency components in the power spectrum of the timing residuals to set a 
95% confidence limit of about 2 x 10~ 8 on flow at a GW frequency of 1/8 yr or 
4 nHz. From Equations [25] and [30l this result is equivalent to a characteristic strain 
at frequency 1/1 yr, Ai yr ss 1.1 x 10 -14 . 

van Haasteren et al. 138 analysed EPTA 1400-MHz data sets for five MSPs with 
spans of 5 - 8 years using a Bayesian analysis to place limits on the GWB amplitude 
as a function of its spectral index a. For a = —2/3, the derived limit at the 95% 
confidence level is Ai yr ss 6 x 10 -15 , about a factor 1.8 better than the Jenet et 
al.® limit. 

NANOGrav multi-band data sets recorded between 2005 and 2010 for 17 MSPs 
were analysed by Demorest et al..® Timing analyses taking into account time- 
varying dispersion delays and frequency-dependent pulse profiles were carried out 
to form sets of post-fit residuals and the corresponding covariance matrices for each 
pulsar. For most MSPs in the sample, no red noise signal was detectable in the 
post-fit residuals. Considering just the pulsar with the smallest post-fit residuals, 
PSR J1713+0747, and taking into account absorption of red noise by the timing 
fit, Demorest et al. obtained a 95% confidence limit for Ai yr of 1.1 x 10 -14 . A sepa¬ 
rate cross-correlation analysis weighted by the expected Hellings & Downs function 
(Equation l22l) across all the pulsars in the sample resulted in a somewhat better 
limit Ai yr «7x 10 -15 , although this limit was dominated by correlations with the 
two best pulsars in the sample, PSRs J1713+0747 and J1909—3744. 

Based on PPTA and earlier Parkes timing observations made in three observing 
bands centred near 700 MHz, 1400 MHz and 3100 MHz respectively with data spans 
of up to 17 years,® together with the PSR B1855+09 archival Arecibo data,® 
Shannon et al. 119 ^ placed a limit of 1.3 x 10 -9 on TLqw at a GW frequency of 2.8 
nHz. The corresponding limit on Ai yr assuming a GW spectral index of —2/3 is 
2.4 x 10 -15 . This analysis, which included dispersion correction for the PPTA data 
sets and was based on the six best PPTA pulsars, used a statistical method similar 
to that of Jenet et al.® but included modelling of the red noise in the timing 
residuals. As shown in Figure l2lil this limit rules out a model in which the growth 
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of SMBH in galaxies is dominated by merger^ a t the 91% confidence level, but 
is consistent with other models for galaxy and SMBH evolution where much of the 
SMBH growth is by accretion. 



Fig. 20. Limits on the relative energy density of the GWB, Qgw a t a GW frequency of 2.8 nHz 
based on the PPTA data sets, together with predictions for PLgw based on several different models 
for the GWB EH The solid and dashed lines that are asymptotic to 1.0 at low TIgw show the 
probability Pr that a GWB signal of energy density f ^gw can exist in the PPTA data sets, based 
on gaussian and non-gaussian GWB statistics respectively. The shaded regio n is ruled out with 95% 
confidence by the PPTA data. Corresponding limits from analysis of EPTA fSSand NANOGrai^S-- 
data sets, scaled to /gw = 2.8 nHz, are also shown. The gaussian curves show the probabilility 
density functions pM for the existence of a GWB with energy density TIgw based on a merger- 
driven model for growth of SMBHs in galaxiespE an empirical synthesis of observational constraints 
on SMBHs in galaxiesU-hii and based on the Millenium dark matter simulations-23 together with 
semi-analytic models for growth of SMBHs in galaxies. (See Ref. 11191 for more detail.) 


As discussed in Section 13.2.21 topological defects in the early Universe are an¬ 
other potential source for the GWB. Figure I2T1 shows limits on the dimensionless 
string tension G/i/c 2 as a function of loop size for various sets of other relevant 
parameters? 11 The middle solid curves are limits based on the current EPTA data 
sets and the lower dashed curves are projections for LEAP data sets that coher¬ 
ently combine data from the EPTA telescopes. The upper dot-dashed line is a limit 
based on LIGO data! 1 The current EPTA results give a conservative upper limit on 
the string tension of 5.3 x 10 -7 . Lower limits can be obtained with more restricted 
assumptions about the string parameters (e.g., Refs 1371 ILi } 1138 ). 
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Fig. 21. Limits on cosmic string tension as a function of the loop size scale parameter a for 
different cutoff node numbers n» and intrinsic spectral slopes q for the current EPTA limit on 
the energy density of the GWB (thick lines) and the projected sensitivity of the LEAP PTA 
(corresponding thin lines). The line with one long dash and three short dashes is an analytic 
approximation which is valid for large loops. The uppermost line is the LIGO limit at / = 1 kHz. 

(Refs mm rq 

3.3.3. Limits on GW emission from individual black-hole binary systems 

For an isolated binary system at a luminosity distance di, the intrinsic GW strain 
amplitude is given by 

, (GM c f/ 3 (nf) 2 / 3 

h ° = — - d^~ (34) 

where M c is the binary chirp mass and / = 2 fb is the GW frequency. The actual 
observed signal depends on the orbital orientation and phase as well as the GW 
polarisation angle. By averaging over these quantities, PTAs can set probabilistic 
limits for the strain amplitude as a function of /, both in a given direction and 
averaged over all directions (see, e.g., Refs [9], d58[). In these analyses, there is as¬ 
sumed to be negligible evolution of / over the data span and only the Earth term 
(cf. Equation 1201) is considered; because of uncertainties in the pulsar distances, the 
pulsar terms cannot be added coherently and just contribute noise. 

Figure (22] shows both sky-averaged upper limits and detection sensitivity for 
continuous-wave GW signals as a function of GW frequency based on the PPTA 
data set and using a frequentist analysis method !- 158 The best limits and sensitivity 
are obtained for GW frequencies around 10 ~ 8 Hz, where the upper limit on ho 
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is about 1.1 x 10 -14 . A similar analysis of the five-year NANOGrav data set by 
Arzoumanian et al. 9 using both frequentist and Bayesian analysis methods gave a 
somewhat higher sky-averaged upper limit of about 5 x 10 -14 at 10 -8 Hz. Because 
of the uneven sky distribution of PTA pulsars, there is quite a strong dependence 
of sensitivity on source direction. This is illustrated in Figure [23] which shows that 
sensitivity is greater toward the greatest concentration of PTA pulsars, roughly 
toward the Galactic Centre. 



Fig. 22. Sky-averaged limits on the intrinsic GW strain amplitude ho as a function of GW 
frequency / based on the PPTA data sets. The lower curves represent the largest GW signal 
(with a false-alarm probability of 1%) that could be present in the real PPTA data (dashed line) 
and a simulated data set (solid line). The upper solid line gives the sensitivity of the PPTA to 
a continuous-wave source, i.e., the minimum signal that could be detected with 95% probability. 
The upper limits and sensitivities are higher at frequencies of 1/1 yr and 1/6 months as these 
frequencies are absorbed by the timing fits for position and parallax respectively. The sloping dot- 
dashed lines are the expected signal levels for a SMBH binary systems with M c = 10 10 Mq and 
distance 400 Mpc (upper line) and M c = 10 9 Mq and distance 30 Mpc (lower line). (Ref. fl58)> 


Figure [22] also shows that we can effectively rule out the existence of SMBH 
binary systems with M c = 10 9 Mq and orbital frequencies around ICC 8 Hz at 
distances closer than 30 Mpc. Similarly, a system with M c = 10 10 Mq at a distance 
of 400 Mpc should be detectable. Unfortunately, as Figure [23] shows, the nearby 
galaxy clusters such as Virgo, Coma and Fornax are all in regions of relatively low 
sensitivity for the PPTA (and other PTAs), so the effective limits for these clusters 
are a factor of a few higher. It is unlikely that such massive binary systems exist 
in these clusters. More generally, limits can also be placed on the SMBH binary 
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Dec=90 



d L (Mpc) Dec=-90° 


15 20 25 50 35 40 45 50 55 60 65 

Fig. 23. Sky distribution of the luminosity distance cLl to which a binary system with chirp mass 
M c = 10 9 Mq radiating at 10 —8 Hz could be detected. The stars indicate the positions of the 20 
PPTA pulsars and the diamonds are potential sources of GW continuous-wave emission. (Ref. [158} 


coalescence rate in the nearby Universe (z < 0.1). Based on the PPTA results, Zhu 
et al.USSI place a 95% confidence limit of 4 x 10 _3 (10 1o Mq/.M c ) 10 / 3 Mpc -3 Gyr -1 
on the coalescence rate. This limit is about two orders of magnitude above current 
estimates of the galaxy merger rate in the local Universe (cf. Ref. [S]) • 

These rate estimates are based on the orientation-averaged and sky-averaged 
amplitudes. It is of course possible that a favourably oriented and located SMBH 
binary system in the late stages of coalescence could exist. On the other hand, the 
estimates are based on circular binary orbits and, as discussed in Sections 13.2.II and 
13.2.31 above, short-period SMBH binaries may have significant eccentricity which 
reduces the GW power at the fundamental frequency / = 2/f, and hence the de¬ 
tectability of such systems. On balance, it seems unlikely that GW from an individ¬ 
ual coalescing SMBH binary system will be detected with the current generation of 
PTAs. 


3 . 4 . Future Prospects 

PTAs have now achieved data spans and ToA precisions that would allow detection 
of the GWB predicted by some models for the evolution of galaxies and the SMBHs 
at their core (see, e.g., Ref. HHJ). Up to now, no detections have been made. While 
this is disappointing from the point of view of GW astrophysics, it is starting to 
have important implications for galaxy and SMBH evolution models and to rule out 
some scenarios. It also implies that PTAs are close to detecting the GWB if current 
predictions for its amplitude are correct. 

Siemens et al. 126 have considered the sensitivity of an idealised PTA to a GWB. 
At low signal levels, when the lowest signal frequencies are below the white noise 
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level, the detection signal-to-noise ratio (S/N) is 


(p) oc Mc-^T? 

G 2, 


(35) 


where M is the number of pulsars in the array, c is the observing cadence (frequency 
of observations), A lyr is the GWB amplitude ('Equation [3(1. a is the rms level of 
the white timing noise, T is the observing data span and /3 is the inverse spectral 
index of the GWB signal in the timing residuals, taken to be 13/3 (Equation I3T1) . 
In the detection regime where the GWB signal exceeds the white noise level, the 
S/N is 



(36) 


Consequently, in the pre-detection regime, the S/N increases rapidly with increased 
observing cadence and data span and decreased timing noise, but has a much weaker 
dependence on these parameters in the strong signal regime. The reason for this is 
that noise from the uncorrelated pulsar term, which is also proportional to Ai yr , 
dominates over the white “receiver” noise, greatly modifying the statistical be¬ 
haviour. Importantly though, in both regimes, the S/N is proportional to M, the 
number of pulsars in the PTA. Figure [Ml illustrates these dependencies for a range 
of plausible future PTAs. 

As Siemens et al4^ point out, if Ai yr ~ 10 -15 , current PTAs are already in 
the “strong signal” regime. This means that increasing the observing data spans 
and cadence or decreasing ToA uncertainties has limited effect on the S/N of a 
potential detection. Increasing the number of pulsars in the PTA is a much more 
cost-effective way to increase detection sensitivity. This fact provides much of the 
motivation to combine data from existing PTAs to form the International Pulsar 
Timing Array ( IPTA ). In the future, the Chinese FAST radio telescoped and the 
SKA 24 will provide a large increase in radiometer sensitivity compared to existing 
instruments. The discussion above shows that this increased sensitivity will be best 
employed in increasing the number of (weaker) pulsars that are timed, rather than 
improving the ToA precision on the stronger existing PTA pulsars. Considerations 
of “jitter noise” in ToAs resulting from shape variations in individual pulsef^lead 
to the same conclusion. 

While direct detection of the GWB would be enormously exciting and signifi¬ 
cant, there is no doubt that direct detection of GW from individual SMBH binary 
systems is potentially much more interesting from an astrophysical perspective. It 
opens up the possibility of identifying a GW source with source or region identified 
through electromagnetic-wave (radio, optical, X-ray or 7 -ray) emission and the ad¬ 
vent of “multi-messenger” astronomy (e.g., Refs l78 | 11141 1221) . For example, many 
active galactic nuclei (AGNs) show evidence of a close binary SMBH at their core. 
Examples of this include X-shaped radio lobes (e.g., Ref. H571) . double-peaked or 
variable emission lines from the core region (e.g., Ref. I124jl . quasi-periodic modula- 
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Fig. 24. Detection S/N for a GWB as a function of PTA data span for four different, but plausible, 
future PTAs. See text for the meaning of the PTA parameters. (Ref. I126D 


tion of core radio (e.g., Ref. (137 1 or X-ray emission (e.g., Ref. ' 1141 [79) . and direct 
imaging of double AGN (e.g., Ref. 1109 ). 

With one exception, existing PTA systems do not have sufficient sensitivity to 
detect these potential GW sources. The exception is the claimed SMBH binary 
system identified by VLBI astrometry of the nearby quasar 3C 66B 133 which was 
effectively ruled out by pulsar timing observations.^ Future PTAs including FAST 
and the SKA will have much increased sensitivity, making searches for other GW 
candidate sources potentially more productive. 

Similarly, although blind searches for continuous-wave GW signals in current 
PTA data sets have a low probability of successful detection, in future this should 
not be the case. The possibility of identification of a source galaxy or AGN then 
depends critically on the accuracy of the position determination for the GW source. 
When only detection of the “Earth term” is considered, the accuracy is at best 
many tens of square degrees (e.g., Ref. 1158,) containing thousands if not millions of 
galaxies. Only a correlation of the GW signal with a modulation of some property 
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of a galaxy core (e.g., intensity or velocity) would establish an identification. 

The situation changes dramatically if the pulsar terms can be added coherently 
with the Earth term (cf. Equation l20l) . Each Earth pulsar system then forms 
an interferometer with baseline d and fringe spacing ~ Xqw /[27rZ(1 + 7 )], where 
d is the pulsar distance, A gw is the GW wavelength and 7 is the direction cosine 
between the pulsar direction and the GW propagation direction . 78 — Positional 
accuracies are roughly the fringe spacing divided by the S/N of the GW detection. 
Since d is typically > 10 3 light-years and A gw is a few light-years, sub-arc-minute 
positional accuracies are possible for the stronger sources. However, to achieve the 
coherent summation, the distance to the pulsars must be known to better than 
Agw/(1 + 7)1 be., about 1 pc unless the GW source is nearly aligned with the 
pulsar. Currently, only one PTA pulsar, PSR J0437—4715, has a distance known to 
this accuracy, measured through VLBI astrometry ,^ 9 but in the SKA era this will 
change. As Boyle and Pen 19 point out, with a high density of PTA pulsars on the 
sky, advantage can be taken of the (1 + 7 ) factor and so pulsars with less precisely 
determined distances located in the general direction of the GW source could be 
useful. In this situation, PTAs would not be confusion-limited and, in principle, 
many individual GW sources could be identified. With a high enough number of 
PTA pulsars (say 1000 or more) it may be possible to localise a binary GW source 
by the quadrupolar pattern of timing residuals in pulsars surrounding the source, 
even if the pulsar distances are poorly determined. 


4. Summary and Conclusions 

Nature has been very kind in providing us with a set of near-perfect celestial clocks, 
many in situations of rapidly varying gravitational accelerations. Not only are these 
celestial clocks, known as pulsars, precise time-keepers, they are also exceedingly 
compact. This enables them to be treated as point masses in theoretical analyses of 
their motion and also permits tests in the regime of strong gravitational fields. These 
qualities result in a very wide range of applications for pulsar time-keeping, most 
importantly, at least in the context of this review, to investigations of relativistic 
gravitation. 

Observations of double-neutron-star systems, wide circular pulsar - white-dwarf 
systems and even isolated pulsars have been used to test the accuracy of gravita¬ 
tional theories. Remarkably, general relativity (GR) is unscathed by all of these tests 
and hence remains the most viable theory of gravitation. Pulsar timing has provided 
the strongest available limits on at least six parameters describing deviations from 
GR. Continued and improved pulsar timing measurements, especially with new and 
highly sensitive radio telescopes such as FAST and the SKA, will both improve on 
these limits and enable new and different tests of relativistic gravity. They may even 
demonstrate a failure of GR to adequately account for the observations, leading to 
new or modified theories of gravitation. 

Continuing and new searches for previously unknown pulsars, especially with 
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FAST and the SKA, will not only increase the number of pulsars that can be used in 
tests of relativistic gravity. They will also turn up new and exciting classes of object 
such as the recently discovered triple system, PSR J0337+1715. Such discoveries 
enrich the investigations of relativistic astrophysics that can be undertaken with 
pulsars. 

One of the outstanding goals of current astronomy and astrophysics is the direct 
detection and study of gravitational waves (GWs). Pulsar timing arrays (PTAs) 
provide a viable mechanism for detection of GWs with frequencies in the nanohertz 
range. They therefore complement other existing or planned instruments such as 
the laser-interferometer systems LIGO and eLISA which are sensitive to GWs at 
frequencies of around 100 Hz and millihertz respectively. The most probable sources 
for GW detection by PTAs are binary super-massive black holes in the cores of dis¬ 
tant galaxies. These produce an unresolved background of GWs that is potentially 
detectable, but there may also be individual binary systems that could be detected 
by PTAs. 

There are currently three major PTAs operating, one each in Europe ( EPTA ), 
North America ( NANOGrav ) and Australia ( PPTA ). Up to now, no GWs have 
been detected by PTAs (or other GW detection systems) so the direct detection 
of GWs remains a goal. However, recent limits on the GW background are placing 
significant constraints on existing models for galaxy mergers over cosmological time 
and the formation and evolution of super-massive black holes in the cores of these 
galaxies. For example, a model in which black-hole growth is dominated by mergers 
is essentially ruled out. 

The sensitivity of PTAs to GWs is a function of several factors including the 
precision of the pulse arrival-time measurements, the data span of the PTA observa¬ 
tions and the cadence or frequency of observations within this data span. However 
the most important single factor is the number of pulsars in the PTA. Of course, 
these pulsars must meet certain timing-precision and period-stability criteria to use¬ 
fully contribute to a PTA. There are two main approaches to increasing the number 
of pulsars. First, existing data sets can be combined to form a single PTA - this 
is the goal of the International Pulsar Timing Array project. Second, searches can 
be undertaken to increase the number of known pulsars suitable for PTA projects. 
With FAST and the SKA it is possible that hundreds of MSPs that are suitable 
for PTA projects will be both discovered and subsequently timed to high precision. 
This will surely lead to the detection of GWs and to detailed investigations of both 
the GWs themselves and the sources that generate them. 
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